simple geometric model

description of all variables

labeled_cartoon.png
\begin{equation} F=\sigma\left(S_{0}+S_{m}+S_{i}\right)+\frac{K}{2} S_{i}\left(\frac{2}{R_{i}}-2C_{i}\right)^{2}+\gamma_{d} S_{d}-\gamma_{m} S_{m}-\gamma_{i} S_{i}+P V_{i} \end{equation}
description value units reference type python variable associated equation LaTeX variable
invagination height   \(nm\)   input (variable) H_i   \(H_{i}\)
membrane tension 0.001 \(pN/nm\) hassinger_design_2017 fixed sigma   \(\sigma\)
total plane radius 2000 \(nm\)   fixed A_t   \(A_{t}\)
final vesicle radius 50 \(nm\) kaksonen_mechanisms_2018 fixed vesicle_radius   \(R_{v}\)
coat bending rigidity 160 \(pN \cdot nm\) dmitrieff_membrane_2015 fixed K   \(K\)
spontaneous curvature 1/50 \(nm^{-1}\) hassinger_design_2017 fixed C_i   \(C_{i}\)
droplet-cytosol surface tension 0.005 \(pN/nm\) jawerth_salt-dependent_2018 fixed gamma_d   \(\gamma_{d}\)
droplet-membrane wetting tension 0.002 \(pN/nm\)   fixed gamma_m   \(\gamma_{m}\)
initial droplet plane radius       fixed      
droplet-invagination surface tension 0.002 \(pN/nm\)   fixed gamma_i   \(\gamma_{i}\)
turgor pressure 0.1 \(pN/nm^2\) ma_mechanics_2019 fixed P   \(P\)
relative invagination height   unitless   input (variable) h rel-h \(h\)
droplet volume   \(nm^3\)   intermediate output V_d   \(V_{d}\)
invagination surface area   \(nm^2\)   intermediate output S_i invagination-area \(S_{i}\)
initial droplet spherical radius 600 \(nm\)   intermediate output drop_radius    
droplet-flat membrane contact area   \(nm^2\)   intermediate output S_m md-area \(S_{m}\)
invagination spherical radius   \(nm\)   intermediate output R_i invagination-radius \(R_{i}\)
droplet-cytosol contact area   \(nm^2\)   intermediate output S_d dc-contact-area \(S_{d}\)
droplet contact angle   rad   intermediate output theta_d droplet-contact-angle \(\theta_{d}\)
invagination contact angle   rad   intermediate output theta_i invagination-angle \(\theta_{i}\)
invagination plane radius   \(nm\)   intermediate output A_i invagination-plane-radius \(A_{i}\)
droplet + invagination plane radius   \(nm\)   intermediate output A_d droplet-plane-radius \(A_{d}\)
outside membrane surface area   \(nm^2\)   intermediate output S_o outside-membrane-area \(S_{o}\)
invagination volume   \(nm^3\)   intermediate output V_i invagination-volume \(V_{i}\)
droplet + invagination volume   \(nm^3\)   intermediate output V_b di-volume \(V_{b}\)
droplet + invagination height   \(nm\)   intermediate output H_d di-height \(H_{d}\)
droplet + invagination spherical radius   \(nm\)   intermediate output R_d di-cap-radius \(R_{d}\)
free energy   \(pN \cdot nm\)   output F f-forces f-shapes \(F\)
membrane tension energy   \(pN \cdot nm\)   output membrane_tension    
bending energy   \(pN \cdot nm\)   output bending_energy    
surface tension energy   \(pN \cdot nm\)   output surface_tension    
wetting energy   \(pN \cdot nm\)   output wetting_energy    
turgor pressure energy   \(pN \cdot nm\)   output turgor_pressure    
membrane curvature   \(nm^{-1}\)   output 2÷R_i    

4.114 \(pN \cdot nm\) = 1 \(k_{B}T\) noauthor_kt_2021

equations

summary

labeled_cartoon.png

arranged by force contributions:

F = membrane tension + bending energy + droplet-cytosol surface tension -
droplet-membrane surface tension - droplet-invagination surface tension + turgor
pressure 
\begin{equation} F=\sigma\left(S_{o}+S_{m}+S_{i}\right)+\frac{K}{2} S_{i}\left(\frac{2}{R_{i}}-2C_{i}\right)^{2}+\gamma_{d} S_{d}-\gamma_{m} S_{m}-\gamma_{i} S_{i}+P V_{i} \end{equation}

rearranged by shapes:

\begin{equation} F= \sigma S_{o} + \left(\sigma - \gamma_{m} \right) S_{m} + \gamma_{d} S_{d} + \frac{K}{2} S_{i} \left(\frac{2}{R_{i}} - 2 C_{i}\right)^{2} + P V_{i} + \left(\sigma - \gamma_{i} \right) S_{i} \end{equation}

broken down to shape functions and constant \(C\):

\begin{equation} F(\Delta S_{o}, \Delta S_{m}, S_{d}, R_{i}, V_{i})= F_{1}(\Delta S_{o}) + F_{2}(\Delta S_{m}) + F_{3}(S_{d}) + F_{4}(R_{i}) + F_{5}(V_{i}) + C \end{equation} \begin{equation} F_{1}(\Delta S_{o}) = \sigma \Delta S_{o} \end{equation} \begin{equation} \Delta S_{o}=\alpha \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} \end{equation} \begin{equation} \label{alpha} \alpha = - \left(\frac{\pi\left(\frac{\gamma_{m}}{\gamma_{d}}+1\right)^{3}} {\left(-\frac{\gamma_{m}}{\gamma_{d}}+1\right)\left(\frac{\gamma_{m}}{\gamma_{d}}+2\right)^{2}}\right)^{1 / 3} \end{equation} \begin{equation} \label{} F_{2}(\Delta S_{m}) = \left(\sigma - \gamma_{m} \right) \Delta S_{m} \end{equation} \begin{equation} \Delta S_{m}= - \alpha \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} + S_{i} h^{2} \end{equation} \begin{equation} \label{} F_{3}(S_{d}) = \gamma_{d} S_{d} \end{equation} \begin{equation} S_{d}=\epsilon \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} \end{equation} \begin{equation} \epsilon = \left(\frac{8 \pi}{\left(-\frac{\gamma_{m}}{\gamma_{d}}+1\right) \left(\frac{\gamma_{m}}{\gamma_{d}}+2\right)^{2}}\right)^{1 / 3} \end{equation} \begin{equation} \label{} \begin{array}{l} F_{4}(R_{i}) = \frac{K}{2} S_{i} \left(\frac{2}{R_{i}} - 2 C_{i}\right)^{2} \end{array} \end{equation} \begin{equation} R_{i} = \frac{\sqrt{S_i}}{2 \sqrt{\pi } h} \end{equation} \begin{equation} \label{} F_{4}(h) = 2 K \left(C_i \sqrt{S_i}-2 \sqrt{\pi } h\right)^{2} \end{equation} \begin{equation} \label{} F_{5}(V_{i}) = P V_{i} \end{equation} \begin{equation} \label{} V_{i}=\eta \left(3 h-2 h^{3}\right) \end{equation} \begin{equation} \label{} \eta=\frac{S_i^{3/2}}{6 \sqrt{\pi }} \end{equation} \begin{equation} \label{} C = \sigma \pi A_t^2+S_i \left(\gamma _m-\gamma _i\right) \end{equation}

all combined into full free energy equation:

\begin{equation} F(h)=\left(\gamma_{m} \alpha+\gamma_{d} \epsilon\right)\left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3}+\left(\sigma-\gamma_{m}\right) S_{i} h^{2} +2 K\left(C_{i} \sqrt{S}_{i}-2 \sqrt{\pi} h\right)^{2}+P \eta \left(3 h-2 h^{3}\right)+C \end{equation}

expanded polynomial:

\begin{equation} F(h)= \left(\gamma_{m} \alpha+\gamma_{d} \epsilon\right) \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} - 2 P \eta h^{3}+\left(\left(\sigma-\gamma_{m}\right) S_{i}+8 K \pi\right) h^{2} +\left(3 P \eta-8 K C_{i} \sqrt{\pi S_{i}}\right) h+B \end{equation} \begin{equation} B=C+2 K C_{i}^{2} S_{i} \end{equation}

individual shape solutions

relative invagination height

h.jpg
\begin{equation} \label{rel-h} h = \frac{H_{i}}{\left(\frac{S_{i}}{\pi}\right)^{1/2}} \end{equation} \begin{equation} H_{i} = h \left(\frac{S_{i}}{\pi}\right)^{1/2} \end{equation}

droplet-membrane contact angle

2020-09-21_19-17-56_screenshot.png
Figure 1: contact angle cartoon contact-angle
\begin{equation} \begin{array}{l} \gamma_{S G}=\gamma_{S L}+\gamma_{L G} \cos \theta \\ \sigma+\gamma_{m}=\sigma+\gamma_{d} \cos \theta_{d} \\ \cos \theta_{d}=\frac{\gamma_{m}}{\gamma_{d}} \\ \sin \theta_{d}=\sqrt{1 - \cos \theta_{d}^{2}} \\ \sin \theta_{d}=\sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}} \end{array} \end{equation} \begin{equation} \label{droplet-contact-angle} \theta_{d}=\arccos \left(\frac{\gamma_{m}}{\gamma_{d}}\right) \end{equation}

outside membrane area

\begin{equation} \label{outside-membrane-area} S_{0}=\pi\left(A_{t}^{2}-A_{d}^{2}\right) \end{equation}

circular plane radii

\begin{equation} \begin{array}{l} \sin \theta=\frac{a}{r} \\ a=r \sin \theta} \\ A_{i}=R_{i} \sin \theta_{i} \end{array} \end{equation} \begin{equation} \label{droplet-plane-radius} A_{d}=R_{d} \sin \theta_{d} \end{equation}

droplet-cytosol contact area

\begin{equation} \label{dc-contact-area} \begin{array}{l} A=2 \pi r^{2}(1-\cos \theta) \\ S_{d}=2 \pi R_{d}^{2}\left(1-\cos \theta_{d}\right) \end{array} \end{equation}

invagination surface area

\begin{equation} \label{invagination-area} S_{i}=4 \pi R_{v}^{2} \end{equation}

invagination volume

\begin{equation} \label{invagination-volume-a-h} \begin{array}{l} V=\frac{1}{6} \pi h\left(3 a^{2}+h^{2}\right) \\ V_{i}=\frac{1}{6} \pi H_{i}\left(3 A_{i}^{2}+H_{i}^{2}\right) \end{array} \end{equation}

alternately:

\begin{equation} \label{invagination-volume} \begin{array}{l} V=\frac{\pi}{3} h^{2} \left(3 r - h \right) \\ V_{i}=\frac{\pi}{3} H_{i}^{2} \left(3 R_{i} - H_{i} \right) \end{array} \end{equation}

invagination plane radius

\begin{equation} \begin{array}{l} A=\pi\left(a^{2}+h^{2}\right) \\ S_{i}=\pi\left(A_{i}^{2}+H_{i}^{2}\right) \\ A_{i}=\sqrt{\frac{S_{i}}{\pi}-H_{i}^{2}} \end{array} \end{equation} \begin{equation} \label{invagination-plane-radius} A_{i}=\sqrt{\frac{S_{i}}{\pi} \left(1 - h^{2}\right)} \end{equation}

invagination radius

\begin{equation} \label{invagination-radius} \begin{aligned} A &=2 \pi r h \\ S_{i} &=2 \pi R_{i} H_{i} \\ R_{i} &=\frac{S_{i}}{2 \pi H_{i}} \end{aligned} \end{equation}

invagination contact angle

\begin{equation} \label{invagination-angle} \begin{array}{l} S_{i}=2 \pi R_{i}^{2}\left(1-\cos \theta_{i}\right) \\ \theta_{i}=\arccos \left(1-\frac{S_{i}}{2 \pi R_{i}^{2}}\right) \end{array} \end{equation}

TODO droplet volume

  • should depend on initial plane radius

    \begin{equation} \label{d-volume} V_{b}=V_{d}+V_{i} \end{equation}

droplet + invagination volume

\begin{equation} \label{di-volume} V_{b}=V_{d}+V_{i} \end{equation}

droplet + invagination spherical cap radius

\begin{equation} \label{di-cap-radius} \begin{array}{l} V=\frac{\pi}{3} r^{3}(2+\cos \theta)(1-\cos \theta)^{2} \\ V_{b}=\frac{\pi}{3} R_{d}^{3}\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2} \\ R_{d}=\left(\frac{3 V_{b}}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{1 / 3} \end{array} \end{equation}

droplet + invagination height

\begin{equation} \label{di-height} \begin{array}{l} A=2\pi r h \\ S_{d}=2\pi R_{d} H_{d} \\ H_{d}=\frac{S_{d}}{2 \pi R_{d}} \end{array} \end{equation}

flat membrane-droplet contact area

\begin{equation} \label{md-area} S_{m}=\pi\left(A_{d}^{2}-A_{i}^{2}\right) \end{equation}

separate force components (original version)

\begin{equation} \label{f-forces} F=\sigma\left(S_{o}+S_{m}+S_{i}\right)+\frac{K}{2} S_{i}\left(\frac{2}{R_{i}}-2C_{i}\right)^{2}+\gamma_{d} S_{d}-\gamma_{m} S_{m}-\gamma_{i} S_{i}+P V_{i} \end{equation}

fully substituted

\begin{equation} \label{free-energy} F=\sigma(\pi(A_{t}^{2}-(((\frac{3 (V_{d}+(\frac{1}{6} \pi \left(\frac{S_{i} h}{\pi}\right)(3 (\sqrt{\frac{(4 \pi R_{v}^{2})}{\pi}-\left(\frac{S_{i} h}{\pi}\right)^{2}})^{2}+\left(\frac{S_{i} h}{\pi}\right)^{2})))}{\pi(2+(\frac{\gamma_{m}}{\gamma_{d}}))(1-(\frac{\gamma_{m}}{\gamma_{d}}))^{2}})^{1 / 3}) (\sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}}))^{2})+(\pi((((\frac{3 (V_{d}+(\frac{1}{6} \pi \left(\frac{S_{i} h}{\pi}\right)(3 (\sqrt{\frac{(4 \pi R_{v}^{2})}{\pi}-\left(\frac{S_{i} h}{\pi}\right)^{2}})^{2}+\left(\frac{S_{i} h}{\pi}\right)^{2})))}{\pi(2+(\frac{\gamma_{m}}{\gamma_{d}}))(1-(\frac{\gamma_{m}}{\gamma_{d}}))^{2}})^{1 / 3}) (\sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}}))^{2}-(\sqrt{\frac{(4 \pi R_{v}^{2})}{\pi}-\left(\frac{S_{i} h}{\pi}\right)^{2}})^{2}))+(4 \pi R_{v}^{2}))+ \frac{K}{2} (4 \pi R_{v}^{2})(\frac{2}{(\frac{(4 \pi R_{v}^{2})}{2 \pi \left(\frac{S_{i} h}{\pi}\right)})}-2C_{i})^{2}+ \gamma_{d} (2 \pi R_{d}^{2}(1-(\frac{\gamma_{m}}{\gamma_{d}})))- \gamma_{m} (\pi((((\frac{3 (V_{d}+(\frac{1}{6} \pi \left(\frac{S_{i} h}{\pi}\right)(3 (\sqrt{\frac{(4 \pi R_{v}^{2})}{\pi}-\left(\frac{S_{i} h}{\pi}\right)^{2}})^{2}+\left(\frac{S_{i} h}{\pi}\right)^{2})))}{\pi(2+(\frac{\gamma_{m}}{\gamma_{d}}))(1-(\frac{\gamma_{m}}{\gamma_{d}}))^{2}})^{1 / 3}) (\sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}}))^{2}-(\sqrt{\frac{(4 \pi R_{v}^{2})}{\pi}-\left(\frac{S_{i} h}{\pi}\right)^{2}})^{2}))- \gamma_{i} (4 \pi R_{v}^{2})+ P (\frac{1}{6} \pi \left(\frac{S_{i} h}{\pi}\right)(3 (\sqrt{\frac{(4 \pi R_{v}^{2})}{\pi}-\left(\frac{S_{i} h}{\pi}\right)^{2}})^{2}+\left(\frac{S_{i} h}{\pi}\right)^{2})) \end{equation}

separate shape force components

\begin{equation} F(\Delta S_{o}, \Delta S_{m}, S_{d}, R_{i}, V_{i})= F_{1}(\Delta S_{o}) + F_{2}(\Delta S_{m}) + F_{3}(S_{d}) + F_{4}(R_{i}) + F_{5}(V_{i}) + C \end{equation} \begin{equation} \label{f-shapes} F= \sigma S_{o} + \left(\sigma - \gamma_{m} \right) S_{m} + \gamma_{d} S_{d} + \frac{K}{2} S_{i} \left(\frac{2}{R_{i}} - 2 C_{i}\right)^{2} + P V_{i} + D \end{equation}

full derivations

note: there may be some residual messyness since I caught a mistake in the sign of \(\gamma_{d}\) in droplet-contact-angle and corrected it after writing all these equations

  • \(F_{1}(\Delta S_{o})\)
    \begin{equation} F_{6}(S_{o}) = \sigma S_{o} \end{equation}

    substituting outside-membrane-area:

    \begin{equation} F_{6}(S_{o}) = \sigma \pi\left(A_{t}^{2}-A_{d}^{2}\right) \end{equation}

    substituting droplet-plane-radius:

    \begin{equation} F_{6}(S_{o}) = \sigma \pi\left(A_{t}^{2}-(R_{d} \sin \theta_{d})^{2}\right) \end{equation}

    substituting di-cap-radius:

    \begin{equation} F_{6}(S_{o}) = \sigma \pi\left(A_{t}^{2}-\left(\left(\frac{3 V_{b}}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{1 / 3} \sin \theta_{d}\right)^{2}\right) \end{equation}

    substituting droplet-contact-angle:

    \begin{equation} F_{6}(S_{o}) = \sigma \pi\left(A_{t}^{2}-\left(\left(\frac{3 V_{b}}{\pi\left(2+\frac{\gamma_{m}}{\gamma_{d}}\right)\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right)^{2}}\right)^{1 / 3} \sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}}\right)^{2}\right) \end{equation}

    substituting di-volume:

    \begin{equation} F_{6}(S_{o}) = \sigma \pi\left(A_{t}^{2}-\left(\left(\frac{3 (V_{d}+V_{i})}{\pi\left(2+\frac{\gamma_{m}}{\gamma_{d}}\right)\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right)^{2}}\right)^{1 / 3} \sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}}\right)^{2}\right) \end{equation}

    substituting invagination-volume:

    \begin{equation} F_{6}(S_{o}) = \sigma \pi\left(A_{t}^{2}-\left(\left(\frac{3 \left(V_{d}+\frac{\pi}{3} H_{i}^{2} \left(3 R_{i} - H_{i} \right)\right)}{\pi\left(2+\frac{\gamma_{m}}{\gamma_{d}}\right)\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right)^{2}}\right)^{1 / 3} \sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}}\right)^{2}\right) \end{equation}

    substituting invagination-radius:

    \begin{equation} F_{6}(S_{o}) = \sigma \pi\left(A_{t}^{2}-\left(\left(\frac{3 \left(V_{d}+\frac{\pi}{3} H_{i}^{2} \left(3 \frac{S_{i}}{2 \pi H_{i}} - H_{i} \right)\right)}{\pi\left(2+\frac{\gamma_{m}}{\gamma_{d}}\right)\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right)^{2}}\right)^{1 / 3} \sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}}\right)^{2}\right) \end{equation}

    “simplified” in mathematica:

    \begin{equation} F_{6}(S_{o}) = \sigma \pi \left(A_t^2+\frac{\left(\gamma _m^2-\gamma _d^2\right) \left(\frac{\gamma _d^2 \left(3 V_d+\pi H_i^2\left(\frac{3 S_i}{2 \pi H_i}-H_i\right)\right)}{\pi \left(\frac{\gamma _m}{\gamma _d}+2\right) \left(\gamma _d-\gamma _m\right)^{2}}\right)^{2/3}}{\gamma _d^2}\right) \end{equation}

    rearranged to isolate constants vs. functions of \(H_{i}\):

    endocytosis phase separation V2 - Page 5.png
    \begin{equation} F_{6}(S_{o}) = \sigma S_{o} \end{equation} \begin{equation} S_{o} = \Delta S_{o} + \pi A_{t}^{2} \end{equation} \begin{equation} F_{6}(S_{o}) = F_{1}(\Delta S_{o}) + \sigma \pi A_{t}^{2} \end{equation} \begin{equation} F_{1}(\Delta S_{o}) = \sigma \Delta S_{o} \end{equation} \begin{equation} \Delta S_{o} = S_{o} - \pi A_{t}^{2} \end{equation} \begin{equation} \Delta S_{o}=\left(\frac{\pi\left(\gamma_{m}^{2}-\gamma_{d}^{2}\right)^{3}}{\left(2 \gamma_{d}+\gamma_{m}\right)^{2}\left(\gamma_{d}-\gamma_{m}\right)^{4}}\right)^{1 / 3}\left(-\pi H_{i}^{3}+\frac{3}{2} S_{i} H_{i}+3 V_{d}\right)^{2 / 3} \end{equation}

    replacing invagination height from rel-h:

    \begin{equation} H_{i} = h \left(\frac{S_{i}}{\pi}\right)^{1/2} \end{equation} \begin{equation} \Delta S_{o}=\left(\frac{\pi\left(\gamma_{m}^{2}-\gamma_{d}^{2}\right)^{3}}{\left(2 \gamma_{d}+\gamma_{m}\right)^{2}\left(\gamma_{d}-\gamma_{m}\right)^{4}}\right)^{1 / 3}\left(-\pi \left(h \left(\frac{S_{i}}{\pi}\right)^{1/2}\right)^{3}+\frac{3}{2} S_{i} h \left(\frac{S_{i}}{\pi}\right)^{1/2}+3 V_{d}\right)^{2 / 3} \end{equation}

    simplified:

    \begin{equation} \Delta S_{o}=\left(\frac{\pi\left(\gamma_{m}^{2}-\gamma_{d}^{2}\right)^{3}}{\left(2 \gamma_{d}+\gamma_{m}\right)^{2}\left(\gamma_{d}-\gamma_{m}\right)^{4}}\right)^{1 / 3}\left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} \end{equation} \begin{equation} \Delta S_{o}=\alpha \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} \end{equation} \begin{equation} \alpha = \left(\frac{\pi\left(\gamma_{m}^{2}-\gamma_{d}^{2}\right)^{3}}{\left(2 \gamma_{d}+\gamma_{m}\right)^{2}\left(\gamma_{d}-\gamma_{m}\right)^{4}}\right)^{1 / 3} \end{equation}
  • \(F_{2}(\Delta S_{m})\)

    gathering from free-energy:

    \begin{equation} \label{} F_{7}(S_{m}) = \left(\sigma - \gamma_{m} \right) S_{m} \end{equation}

    from md-area:

    \begin{equation} \label{md-area} S_{m}=\pi\left(A_{d}^{2}-A_{i}^{2}\right) \end{equation}

    substitute droplet-plane-radius:

    \begin{equation} A_{d}=R_{d} \sin \theta_{d} \end{equation}

    and invagination-plane-radius:

    \begin{equation} A_{i}=\sqrt{\frac{S_{i}}{\pi}-H_{i}^{2}} \end{equation}

    to get:

    \begin{equation} S_{m}=\pi\left(\left(R_{d} \sin \theta_{d}\right)^{2}-\frac{S_{i}}{\pi}+H_{i}^{2}\right) \end{equation}

    rearranged:

    \begin{equation} S_{m}=\pi\left(\left(R_{d} \sin \theta_{d}\right)^{2}+H_{i}^{2}\right) - S_{i} \end{equation}

    substitute di-cap-radius:

    \begin{equation} R_{d}=\left(\frac{3 V_{b}}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{1 / 3} \end{equation}

    to get:

    \begin{equation} S_{m}=\pi\left(\left(\left(\frac{3 V_{b}}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{1 / 3} \sin \theta_{d}\right)^{2}+H_{i}^{2}\right) - S_{i} \end{equation}

    substitute di-volume:

    \begin{equation} V_{b}=V_{d}+V_{i} \end{equation}

    to get:

    \begin{equation} S_{m}=\pi\left(\left(\left(\frac{3 \left(V_{d}+V_{i} \right)}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{1 / 3} \sin \theta_{d}\right)^{2}+H_{i}^{2}\right) - S_{i} \end{equation}

    substitute invagination-volume:

    \begin{equation} V_{i}=\frac{\pi}{3} H_{i}^{2} \left(3 R_{i} - H_{i} \right) \end{equation}

    to get:

    \begin{equation} S_{m}=\pi\left(\left(\left(\frac{3 \left(V_{d}+\frac{\pi}{3} H_{i}^{2} \left(3 R_{i} - H_{i} \right) \right)}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{1 / 3} \sin \theta_{d}\right)^{2}+H_{i}^{2}\right) - S_{i} \end{equation}

    substitute invagination-radius:

    \begin{equation} R_{i} =\frac{S_{i}}{2 \pi H_{i}} \end{equation}

    to get:

    \begin{equation} S_{m}=\pi\left(\left(\left(\frac{3 \left(V_{d}+\frac{\pi}{3} H_{i}^{2} \left(3 \frac{S_{i}}{2 \pi H_{i}} - H_{i} \right) \right)}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{1 / 3} \sin \theta_{d}\right)^{2}+H_{i}^{2}\right) - S_{i} \end{equation}

    substitute droplet-contact-angle:

    \begin{equation} \begin{array}{l} \cos \theta_{d}=\frac{\gamma_{m}}{\gamma_{d}} \\ \sin \theta_{d}=\sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}} \end{array} \end{equation}

    to get:

    \begin{equation} S_{m}=\pi\left(\left(\left(\frac{3 \left(V_{d}+\frac{\pi}{3} H_{i}^{2} \left(3 \frac{S_{i}}{2 \pi H_{i}} - H_{i} \right) \right)}{\pi\left(2+\frac{\gamma_{m}}{\gamma_{d}}\right)\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right)^{2}}\right)^{1 / 3} \sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}}\right)^{2}+H_{i}^{2}\right) - S_{i} \end{equation}

    simplified in mathematica:

    \begin{equation} S_{m}=\pi \left(\frac{\left(\gamma _d^2-\gamma _m^2\right) \left(\frac{\gamma _d^2 \left(3 V_d+\pi H_i^2\left(\frac{3 S_i}{2 \pi H_i}-H_i\right)\right)}{\pi \left(2-\frac{\gamma _m}{\gamma _d}\right) \left(\gamma _d+\gamma _m\right)^{2}}\right)^{2/3}}{\gamma _d^2}+H_i^2\right)-S_i \end{equation}

    manually rearranged:

    endocytosis phase separation V2 - Page 6.png
    \begin{equation} \label{} F_{7}(S_{m}) = \left(\sigma - \gamma_{m} \right) S_{m} \end{equation} \begin{equation} S_{m}= \Delta S_{m} - S_{i} \end{equation} \begin{equation} \label{} F_{7}(S_{m}) = F_{2}(\Delta S_{m}) - \left(\sigma - \gamma_{m} \right) S_{i} \end{equation} \begin{equation} \label{} F_{2}(\Delta S_{m}) = \left(\sigma - \gamma_{m} \right) \Delta S_{m} \end{equation} \begin{equation} \Delta S_{m}=\left(\frac{\pi\left(\gamma_{d}^{2}-\gamma_{m}^{2}\right)^{3}}{\left(2 \gamma_{d}+\gamma_{m}\right)^{2}\left(\gamma_{d}-\gamma_{m}\right)^{4}}\right)^{1 / 3}\left(-\pi H_{i}^{3}+\frac{3}{2} S_{i} H_{i}+3 V_{d}\right)^{2 / 3}+\pi H_{i}^{2} \end{equation}

    replacing invagination height from rel-h:

    \begin{equation} H_{i} = h \left(\frac{S_{i}}{\pi}\right)^{1/2} \end{equation} \begin{equation} \Delta S_{m}=\left(\frac{\pi\left(\gamma_{d}^{2}-\gamma_{m}^{2}\right)^{3}}{\left(2 \gamma_{d}+\gamma_{m}\right)^{2}\left(\gamma_{d}-\gamma_{m}\right)^{4}}\right)^{1 / 3}\left(-\pi \left(h \left(\frac{S_{i}}{\pi}\right)^{1/2}\right)^{3}+\frac{3}{2} S_{i} \left(h \left(\frac{S_{i}}{\pi}\right)^{1/2}\right)+3 V_{d}\right)^{2 / 3}+ \pi \left(h \left(\frac{S_{i}}{\pi}\right)^{1/2}\right)^{2} \end{equation}

    simplified:

    \begin{equation} \Delta S_{m}=\left(\frac{\pi\left(\gamma_{d}^{2}-\gamma_{m}^{2}\right)^{3}}{\left(2 \gamma_{d}+\gamma_{m}\right)^{2}\left(\gamma_{d}-\gamma_{m}\right)^{4}}\right)^{1 / 3} \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} + S_{i} h^{2} \end{equation} \begin{equation} \Delta S_{m}= \beta \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} + S_{i} h^{2} \end{equation} \begin{equation} \beta = \left(\frac{\pi\left(\gamma_{d}^{2}-\gamma_{m}^{2}\right)^{3}}{\left(2 \gamma_{d}+\gamma_{m}\right)^{2}\left(\gamma_{d}-\gamma_{m}\right)^{4}}\right)^{1 / 3} \end{equation}
  • \(F_{3}(S_{d})\)

    gathering from free-energy:

    \begin{equation} \label{} F_{3}(S_{d}) = \gamma_{d} S_{d} \end{equation}

    from dc-contact-area:

    \begin{equation} S_{d}=2 \pi R_{d}^{2}\left(1-\cos \theta_{d}\right) \end{equation}

    substitute di-cap-radius

    \begin{equation} R_{d}=\left(\frac{3 V_{b}}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{1 / 3} \end{equation}

    to get:

    \begin{equation} S_{d}=2 \pi \left(\frac{3 V_{b}}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{2 / 3}\left(1-\cos \theta_{d}\right) \end{equation}

    substitute di-volume

    \begin{equation} V_{b}=V_{d}+V_{i} \end{equation}

    to get:

    \begin{equation} S_{d}=2 \pi \left(\frac{3 \left(V_{d}+V_{i}\right)}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{2 / 3}\left(1-\cos \theta_{d}\right) \end{equation}

    substitute invagination-volume

    \begin{equation} V_{i}=\frac{\pi}{3} H_{i}^{2} \left(3 R_{i} - H_{i} \right) \end{equation}

    to get:

    \begin{equation} S_{d}=2 \pi \left(\frac{3 \left(V_{d}+\frac{\pi}{3} H_{i}^{2} \left(3 R_{i} - H_{i} \right)\right)}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{2 / 3}\left(1-\cos \theta_{d}\right) \end{equation}

    substitute droplet-contact-angle

    \begin{equation} \cos \theta_{d}=\frac{\gamma_{m}}{\gamma_{d}} \end{equation}

    to get:

    \begin{equation} S_{d}=2 \pi \left(\frac{3 \left(V_{d}+\frac{\pi}{3} H_{i}^{2} \left(3 R_{i} - H_{i} \right)\right)}{\pi\left(2+\frac{\gamma_{m}}{\gamma_{d}}\right)\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right)^{2}}\right)^{2 / 3}\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right) \end{equation}

    substitute invagination-radius

    \begin{equation} R_{i} &=\frac{S_{i}}{2 \pi H_{i}} \end{equation}

    to get:

    \begin{equation} S_{d}=2 \pi \left(\frac{3 \left(V_{d}+\frac{\pi}{3} H_{i}^{2} \left(3 \frac{S_{i}}{2 \pi H_{i}} - H_{i} \right)\right)}{\pi\left(2+\frac{\gamma_{m}}{\gamma_{d}}\right)\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right)^{2}}\right)^{2 / 3}\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right) \end{equation}

    simplified in mathematica:

    \begin{equation} S_{d}= \frac{2 \pi \left(\gamma _d-\gamma _m\right) \left(\frac{\gamma _d^2 \left(3 V_d+\pi H_i^2\left(\frac{3 S_i}{2 \pi H_i}-H_i\right)\right)}{\pi \left(\frac{\gamma _m}{\gamma _d}+2\right) \left(\gamma _d-\gamma _m\right)^{2}}\right)^{2/3}}{\gamma _d} \end{equation}

    manually rearranged:

    20210104_sd.png
    \begin{equation} S_{d}=\left(\frac{8 \pi}{\left(2+\frac{\gamma_{m}}{\gamma_{d}}\right)^{2}\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right)}\right)^{1 / 3}\left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} \end{equation} \begin{equation} S_{d}=\epsilon \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} \end{equation} \begin{equation} \epsilon = \left(\frac{8 \pi}{\left(2+\frac{\gamma_{m}}{\gamma_{d}}\right)^{2}\left(1-\frac{\gamma_{d}}{\gamma_{m}}\right)}\right)^{1 / 3} \end{equation}
  • \(F_{4}(R_{i})\)

    gathering from f-forces:

    \begin{equation} \label{} \begin{array}{l} F_{4}(R_{i}) = \frac{K}{2} S_{i} \left(\frac{2}{R_{i}} - 2 C_{i}\right)^{2} \end{array} \end{equation}

    from invagination-radius

    \begin{equation} R_{i} &=\frac{S_{i}}{2 \pi H_{i}} \end{equation}

    replacing invagination height from rel-h:

    \begin{equation} H_{i} = h \left(\frac{S_{i}}{\pi}\right)^{1/2} \end{equation} \begin{equation} R_{i} = \frac{S_{i}}{2 \pi h \left(\frac{S_{i}}{\pi}\right)^{1/2}} \end{equation}

    simplified:

    \begin{equation} R_{i} = \frac{\sqrt{S_i}}{2 \sqrt{\pi } h} \end{equation}
    • \(F_{4}(h)\)

      substitute invagination-radius

      \begin{equation} R_{i} &=\frac{S_{i}}{2 \pi H_{i}} \end{equation}

      to get:

      \begin{equation} \label{} \begin{array}{l} F_{4}(H_{i}) = \frac{K}{2} S_{i} \left(\frac{2}{\frac{S_{i}}{2 \pi H_{i}}} - 2 C_{i}\right)^{2} \end{array} \end{equation}

      simplified:

      \begin{equation} \label{} F_{4}(H_{i}) = \frac{2 K}{S_i} \left(C_i S_i-2 \pi H_i\right)^{2} \end{equation}

      replacing invagination height from rel-h:

      \begin{equation} H_{i} = h \left(\frac{S_{i}}{\pi}\right)^{1/2} \end{equation} \begin{equation} \label{} F_{4}(h) = \frac{2 K}{S_i} \left(C_i S_i-2 \pi \left(h \left(\frac{S_{i}}{\pi}\right)^{1/2}\right)\right)^{2} \end{equation}

      simplified:

      \begin{equation} \label{} F_{4}(h) = 2 K \left(C_i \sqrt{S_i}-2 \sqrt{\pi } h\right)^2 \end{equation}
  • \(F_{5}(V_{i})\)

    gathering from f-forces:

    \begin{equation} \label{} F_{5}(V_{i}) = P V_{i} \end{equation}

    from invagination-volume-a-h:

    \begin{equation} \label{} \begin{array}{l} V_{i}=\frac{1}{6} \pi H_{i}\left(3 A_{i}^{2}+H_{i}^{2}\right) \end{array} \end{equation}

    substitute invagination-plane-radius

    \begin{equation} \label{} \begin{array}{l} A_{i}=\sqrt{\frac{S_{i}}{\pi}-H_{i}^{2}} \end{array} \end{equation}

    to get:

    \begin{equation} \label{} \begin{array}{l} V_{i}=\frac{1}{6} \pi H_{i}\left(3 \left(\frac{S_{i}}{\pi}-H_{i}^{2}\right)+H_{i}^{2}\right) \end{array} \end{equation}

    simplified:

    \begin{equation} \label{} \begin{array}{l} V_{i}=\frac{1}{6} \pi H_{i}\left(\frac{3 S_i}{\pi }-2 H_{i}^2\right) \end{array} \end{equation}

    replacing invagination height from rel-h:

    \begin{equation} H_{i} = h \left(\frac{S_{i}}{\pi}\right)^{1/2} \end{equation} \begin{equation} \label{} V_{i}=\frac{1}{6} \pi h \left(\frac{S_{i}}{\pi}\right)^{1/2}\left(\frac{3 S_i}{\pi }-2 \left(h \left(\frac{S_{i}}{\pi}\right)^{1/2}\right)^2\right) \end{equation}

    simplified:

    \begin{equation} \label{} V_{i}= \frac{h \left(3-2 h^2\right) S_i^{3/2}}{6 \sqrt{\pi }} \end{equation} \begin{equation} \label{} V_{i}=\eta \left(3 h-2 h^{3}\right) \end{equation} \begin{equation} \label{} \eta=\frac{S_i^{3/2}}{6 \sqrt{\pi }} \end{equation}
  • \(D\): everything else that doesn’t change with \(H_{i}\)
    \begin{equation} \label{D} D = \sigma S_{i} - \gamma_{i} S_{i} \end{equation}
  • \(C\): additional constants grouped out from \(\Delta\) terms
    \begin{equation} \label{} C = D + \sigma \pi A_{t}^{2} - \left(\sigma - \gamma_{m} \right) S_{i} \end{equation} \begin{equation} \label{} D = \sigma S_{i} - \gamma_{i} S_{i} \end{equation}

    simplified:

    \begin{equation} \label{} C = \sigma \pi A_t^2+S_i \left(\gamma _m-\gamma _i\right) \end{equation}

DONE substituted after isolating \(h\)

  • could anything be factored out more first?
\begin{equation} \label{} F(h) = \sigma \alpha \xi + \left(\sigma - \gamma_{m} \right) \left( \beta \xi + S_{i} h^{2}\right) + \gamma_{d} \epsilon \xi + 2 K \left(C_i \sqrt{S_i}-2 \sqrt{\pi } h\right)^{2} + P \eta \left(3 h-2 h^{3}\right) + C \end{equation} \begin{equation} \label{} \xi = \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} \end{equation}
20201203_fh (1).png
\begin{equation} F(h)=\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \epsilon\right)\left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3}+\left(\sigma-\gamma_{m}\right) S_{i} h^{2} +2 K\left(C_{i} \sqrt{S}_{i}-2 \sqrt{\pi} h\right)^{2}+P \eta \left(3 h-2 h^{3}\right)+C \end{equation} \begin{equation} F(h)= \left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \epsilon\right) \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} - 2 P \eta h^{3}+\left(\left(\sigma-\gamma_{m}\right) S_{i}+8 K \pi\right) h^{2} +\left(3 P \eta-8 C_{i} \sqrt{\pi S_{i}}\right) h+B \end{equation} \begin{equation} B=C+2 K C_{i}^{2} S_{i} \end{equation}

testing sign of alpha terms

import numpy as np
import math
def cube_root(i):
    return math.pow(abs(i),float(1)/3) * (1,-1)[i<0]

def alpha_so(gamma_m, gamma_d):
    alpha = cube_root((np.pi*(gamma_m**2 - gamma_d**2)**3)/((2*gamma_d+gamma_m)**2 * (gamma_d-gamma_m)**4))
    return alpha

for gamma_m in np.linspace(0, 200, 9):
    for gamma_d in np.linspace(0, 200, 9):
        print('gamma_m',gamma_m,'gamma_d',gamma_d,alpha_so(gamma_m, gamma_d))

implementation and analysis [36/42]

analytical derivations

thermodynamic spontaneity (global stability) [5/5]

  • I’m defining “thermodynamic spontaneity” as whether F(h=0) > F(h=1)
  • analyzed by plotting the phase boundary line across 2 parameters when F(h=0) == F(h=1)
  • DONE setting up equality
    20210104_thermodynamic_spontaneity_kfix.png
    \begin{equation} \left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \in\right)\left(3 V_{d}\right)^{2 / 3}+B=\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \in\right)\left(3 V_d+\frac{S_i^{3/2}}{2 \sqrt{\pi }}\right)^{2/3} + \left(\sigma-\gamma_{m}\right) S_{i}+8 K \pi+ P \eta-8 K C_{i} \sqrt{\pi S_{i}}+B \end{equation}
  • DONE solving for \(K\)
    20210104_thermodynamic_k.png
    \begin{equation} \label{therm_spont_k} K = \frac{\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \in\right)\left(\left(3 V_{d}\right)^{2 / 3}-\left(\frac{S_{i}^{3 / 2}}{2 \sqrt{\pi}}+3 V_{d}\right)^{2 / 3}\right)-\left(\sigma-\gamma_{m}\right) S_{i}- P \eta}{8 \pi-8 C_{i} \sqrt{\pi S_{i}}} \end{equation}
  • DONE solving for \(P\)
    20210104_thermodynamic_p.png
    \begin{equation} \label{therm_spont_p} P=\frac{\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \in\right)\left(\left(3 V_{d}\right)^{2 / 3}-\left(\frac{S_{i}^{3 / 2}}{2 \sqrt{\pi}}+3 V_{d}\right)^{2 / 3}\right)-\left(\sigma-\gamma_{m}\right) S_{i}-8 K \left(\pi - C_{i} \sqrt{\pi S_{i}}\right)}{\eta} \end{equation}
  • DONE solving for \(\sigma\)
    20210104_thermodynamic_sigma.png
    \begin{equation} \label{therm_spont_sigma} \sigma=\frac{\frac{-\gamma_{m} S_{i}+8 K \left(\pi - C_{i} \sqrt{\pi S_{i}}\right) + P \eta}{\left(3 V_{d}\right)^{2 / 3}-\left(\frac{S_{i}^{3 / 2}}{2 \sqrt{\pi}}+3 V_{d}\right)^{2 / 3}}-\gamma_{d} \epsilon+\gamma_{m} \beta}{\alpha+\beta-\frac{S_{i}}{\left(3 V_{d}\right)^{2 / 3}-\left(\frac{S_{i}^{3 / 2}}{2 \sqrt{\pi}}+3 V_{d}\right)^{2 / 3}}} \end{equation}
  • DONE solving for \(C_{i}\)
    20210104_thermodynamic_ci.png
    \begin{equation} \label{therm_spont_ci} C_{i} = \frac{\pi-\frac{\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \in\right)\left(\left(3 V_{d}\right)^{2 / 3}-\left(\frac{S_{i}^{3 / 2}}{2 \sqrt{\pi}}+3 V_{d}\right)^{2 / 3}\right)-\left(\sigma-\gamma_{m}\right) S_{i}-P_{\eta}}{8 K}}{\sqrt{\pi S_{i}}} \end{equation}
  • any others?

kinetic spontaneity (local stability) [11/11]

  • I’m defining “kinetic spontaneity” as whether F’(h=0) is negative or positive
  • DONE solving derivative
    20210105_f_derivative.png
    \begin{equation} \label{f-prime} F^{\prime}(h)=\frac{\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \epsilon\right)\left(1-2 h^{2}\right) S_{i}^{3 / 2}}{\sqrt{\pi}\left(\frac{\left(\frac{3 h}{2}-h^{3}\right) S_{i}^{3 / 2}}{\sqrt{\pi}}+3 V_{d}\right)^{1 / 3}}-6 P \eta h^{2}+2\left(\left(\sigma-\gamma_{m}\right) S_{i}+8 K \pi\right) h+3 P \eta-8 K C_{i} \sqrt{\pi S_{i}} \end{equation}
  • DONE solving for h=0
    \begin{equation} \label{f-prime-alpha-beta} F^{\prime}(0)=\frac{\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \epsilon\right) S_{i}^{3 / 2}}{\sqrt{\pi} \left(3 V_{d}\right)^{1 / 3}}+3 P \eta-8 K C_{i} \sqrt{\pi S_{i}} \end{equation} \begin{equation} \label{f-prime} F^{\prime}(0)=\frac{\left(\gamma_{m} \alpha + \gamma_{d} \epsilon \right) S_{i}^{3 / 2}}{\sqrt{\pi} \left(3 V_{d}\right)^{1 / 3}}+3 P \eta-8 K C_{i} \sqrt{\pi S_{i}} \end{equation}
  • DONE sign of droplet prefactor at h=0
    20210217_prefactor_sign.png

    the sign of the term concerning the droplet must always be positive.

  • DONE solving for \(K\)
    20210106_kinetic_k.png
    \begin{equation} \label{kinetic_spont_k} K=\frac{\frac{\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \epsilon\right) S_{i}^{3 / 2}}{\sqrt{\pi}\left(3 V_{d}\right)^{1 / 3}}+3 P \eta}{8 C_{i} \sqrt{\pi S_{i}}} \end{equation}
  • DONE solving for \(P\)
    20210106_kinetic_p.png
    \begin{equation} \label{kinetic_spont_p} P=\frac{\frac{-\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \epsilon\right) S_{i}^{3 / 2}}{\sqrt{\pi}\left(3 V_{d}\right)^{1 / 3}}+8 K C_{i} \sqrt{\pi S_{i}}}{3 \eta} \end{equation}
  • DONE solving for \(\sigma\)
    20210106_kinetic_sigma.png
    \begin{equation} \label{kinetic_spont_sigma} \sigma=\frac{\frac{\sqrt{\pi}\left(3 V_{d}\right)^{1 / 3}}{S_{i}^{3 / 2}}\left(8 K C_{i} \sqrt{\pi S_{i}}-3 P \eta \right)-\gamma_{d} \epsilon+\gamma_{m} \beta}{(\alpha+\beta)} \end{equation}
  • DONE solving for \(C_{i}\)
    20210106_derivs_kinetic_ci.png
    \begin{equation} \label{kinetic_spont_ci} C_{i} = \frac{\frac{\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \epsilon\right) S_{i}^{3 / 2}} {\sqrt{\pi}\left(3 V_{d}\right)^{1 / 3}}+3 P \eta}{8 K \sqrt{\pi S_{i}}} \end{equation}
  • CANCELLED solving for \(\gamma_{m}\)
    20210106_derivs_kinetic_gammam.png
    \begin{equation} \label{kinetic_spont_gammam} \gamma_{m}=\frac{\frac{-\sqrt{\pi}\left(3 V_{d}\right)^{1 / 3}}{S_{i}^{3 / 2}}\left(8 K C_{i} \sqrt{\pi S_{i}}-3 P \eta \right)+\gamma_{d} \epsilon+\sigma \alpha+\sigma \beta}{\beta} \end{equation}
  • CANCELLED solving for \(\gamma_{d}\)
    20210106_derivs_kinetic_gammad.png
    \begin{equation} \label{kinetic_spont_gammad} \gamma_{d}=\frac{\frac{\sqrt{\pi}\left(3 V_{d}\right)^{1 / 3}}{S_{i}^{3 / 2}}\left(8 K C_{i} \sqrt{\pi S_{i}}-3 P \eta \right)-\sigma \alpha-\left(\sigma-\gamma_{m}\right) \beta}{\epsilon} \end{equation}
  • CANCELLED solving for \(S_{i}\)
  • DONE solving for \(V_{d}\)
    20210106_derivs_kinetic_vd.png
    \begin{equation} \label{kinetic_spont_vd} V_{d}=\frac{\left(\frac{\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \epsilon \right) S_{i}^{3 / 2}}{\left(8 K C_{i} \sqrt{\pi S_{i}}-3 P \eta \right) \sqrt{\pi}}\right)^{3}}{3} \end{equation}

python code [20/26]

configuration [12/13]

  • import libraries, settings
    import napari
    import os
    from decimal import Decimal
    import math
    import numpy as np
    import pandas as pd
    import matplotlib.pyplot as plt  # plotting
    import matplotlib.animation as animation
    import mpl_toolkits.mplot3d.axes3d as axes3d
    import textwrap
    import seaborn as sns
    import scipy
    #sns.set(style='ticks')
    #sns.set_palette(sns.cubehelix_palette(10, start=.1, rot=-.75, reverse=True))
    sns.set_palette('hls',8)
    #sns.set_palette(sns.diverging_palette(250, 30, l=65, center="dark"))
    # sns.set_palette("coolwarm", 6)
    # from cycler import cycler
    SMALL_SIZE = 12
    MEDIUM_SIZE = 16
    BIGGER_SIZE = 20
    
    
    plt.rc('axes', titlesize=MEDIUM_SIZE)     # fontsize of the axes title
    plt.rc('axes', labelsize=MEDIUM_SIZE)    # fontsize of the x and y labels
    plt.rc('xtick', labelsize=SMALL_SIZE)    # fontsize of the tick labels
    plt.rc('ytick', labelsize=SMALL_SIZE)    # fontsize of the tick labels
    plt.rc('legend', fontsize=SMALL_SIZE)    # legend fontsize
    plt.rc('figure', titlesize=BIGGER_SIZE)  # fontsize of the figure title
    # default_cycler = (cycler(color=['indigo', 'darkgreen', 'purple', 'green', 'magenta', 'lime']))
    # plt.rc('axes', prop_cycle=default_cycler)
    
    #if kernel was started in notebook directory
    #fig_dir = 'figures/'
    #assuming kernel was started in home directory
    fig_dir = 'google_drive/grad_school/db_lab/code/analysis/llps_endocytosis_modeling/figures/'
    animation_dir = 'google_drive/grad_school/db_lab/code/analysis/llps_endocytosis_modeling/animation/'
    
    %matplotlib inline
    %gui qt
    import IPython
    from tabulate import tabulate
    
    class OrgFormatter(IPython.core.formatters.BaseFormatter):
        def __call__(self, obj):
            try:
                return tabulate(obj, headers='keys',
                                tablefmt='orgtbl', showindex='always')
            except:
                return None
    
    ip = get_ipython()
    ip.display_formatter.formatters['text/org'] = OrgFormatter()
    
  • set default fixed parameter values
    • assigned here
      pi = np.pi
      vesicle_radius = 50. #nm
      vesicle_diameter = 2*vesicle_radius #nm
      drop_radius = 1500. #nm
      S_i = 4.*pi*vesicle_radius**2 #nm^2
      #H_i_range = np.linspace(0.001*vesicle_radius,2*vesicle_radius,1000) #nm
      h_range = np.linspace(0.001,1,1000) #dimensionless 0-1
      H_i_range = h_range*(S_i/pi)**(1/2)
      sigma = 150. #arbitrary, N/nm^2
      A_t = 1000 #nm, arbitrary
      K = 100. #arbitrary, N/nm
      C_i = 1/vesicle_radius #nm^-1
      gamma_d = 200. #arbitrary N/nm^2
      gamma_m = 175. #arbitrary, N/nm^2
      gamma_i = gamma_m
      V_d = 4.*pi*drop_radius**2 #nm^3
      P = 1000/(10**9) #N/nm^3
      
      params = {'vesicle_radius' : vesicle_radius,
          'vesicle_diameter' : vesicle_diameter,
          'drop_radius' : drop_radius,
          'S_i' : S_i,
          'h_range' : h_range,
          'H_i_range' : H_i_range,
          'sigma' : sigma,
          'A_t' : A_t,
          'K' : K,
          'C_i' : C_i,
          'gamma_d' : gamma_d,
          'gamma_m' : gamma_m,
          'gamma_i' : gamma_i,
          'V_d' : V_d,
          'P':P}
      
    • assigned from org table

      established in variables

      variable_df = pd.DataFrame(variable_list[1:], columns=variable_list[0])
      variable_df = variable_df.set_index('python variable')
      param_values = variable_df['value']
      variable_df
      
      
      vesicle_radius = param_values['vesicle_radius']
      vesicle_diameter = 2*vesicle_radius #nm
      drop_radius = param_values['drop_radius'] #nm
      S_i = 4.*pi*vesicle_radius**2 #nm^2
      #H_i_range = np.linspace(0.001*vesicle_radius,2*vesicle_radius,1000) #nm
      h_range = np.linspace(0.001,1,1000) #dimensionless 0-1
      H_i_range = h_range*(S_i/pi)**(1/2)
      sigma = param_values['sigma'] #arbitrary, N/nm^2
      A_t = param_values['A_t'] #nm, arbitrary
      K = param_values['K'] #arbitrary, N/nm
      C_i = 1/vesicle_radius #nm^-1
      gamma_d = param_values['gamma_d'] #arbitrary N/nm^2
      gamma_m = param_values['gamma_m'] #arbitrary, N/nm^2
      gamma_i = gamma_m
      #V_d = 4.*pi*drop_radius**2 #nm^3
      V_d = 4/3*pi*drop_radius**3 #nm^3
      P = param_values['P'] #N/nm^2
      
      params = {'vesicle_radius' : vesicle_radius,
          'vesicle_diameter' : vesicle_diameter,
          'drop_radius' : drop_radius,
          'S_i' : S_i,
          'h_range' : h_range,
          'H_i_range' : H_i_range,
          'sigma' : sigma,
          'A_t' : A_t,
          'K' : K,
          'C_i' : C_i,
          'gamma_d' : gamma_d,
          'gamma_m' : gamma_m,
          'gamma_i' : gamma_i,
          'V_d' : V_d,
          'P':P}
      
    • complemented from continuum model

      this is a bit adjusted so that there is a barrier due to the droplet

      N=10
      coat_radius = 50
      #coat_area = ((3*N)/(3*N+2))*4*np.pi*coat_radius**2
      coat_area = 4*np.pi*coat_radius**2
      total_radius = np.sqrt(coat_area/np.pi)*3
      Vd = 1/2*np.pi*200**3
      
      model_params = {
          #'time':40000,
          'time':1000000,
          'dt':1,
          'print_interval':10000,
          'write_interval':1000,
          'shape':'final',
          'P':0.0001,
          'sigma':0.01,
          'gamma_1':0,
          'gamma_2':0.05,
          'gamma_3':0.1,
          'gamma_d':0.2,
          'h1':0.00001,
          'h2':0.00001,
          'h3':0.00001,
          'H1':1/coat_radius,
          'K1':50,
          'H3':0,
          'K3':1,
          'polar_type':'integral',
          'G':1,
          'cell_wall':False,
          'N':N,
          'A_t':total_radius,
          'A_c':coat_area,
          'Vd':Vd,
          'emap':False,
          'tip_pull':0
      }
      
      V_d_range = np.linspace(5E6,5E7,100)
      
  • define composite variables

    This is just to check what they are, I re-define them whenever I solve for F

    alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
    gamma_ratio = gamma_m/gamma_d
    alpha_refactored = np.cbrt((pi*(gamma_ratio+1)**3)/((gamma_ratio-1)*(gamma_ratio+2)**2))
    beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
    #epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
    epsilon = np.cbrt((8*pi)/((2+gamma_ratio)**2*(1-gamma_ratio)))
    eta = (S_i**(3/2))/(6*pi**(1/2))
    C = sigma*pi*A_t**2+S_i*(gamma_m-gamma_i)
    
    alpha, alpha_refactored, beta, epsilon, eta, C
    
  • custom functions [12/13]
    • cube root
      def cube_root(i):
          return math.pow(abs(i),float(1)/3) * (1,-1)[i<0]
      
    • model functions
      • pre-analytical version
        • [ ] rewrite each function variable as dictionary key (outpus['theta_d'] instead of theta_d for the first one) and just return outputs at end
          • should eliminate some redundancy
        def llps_model(H_i=H_i_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                      C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d,
                      P=P):
            #H_i = h*(S_i/pi)**(1/2)
            # might need to be -gamma_m/gamma_d
            theta_d = np.arccos(gamma_m/gamma_d) #equation 8
            #to make a full array out of the constant theta_d value
            #theta_d = np.full_like(H_i, theta_d)
            A_i = ((S_i/pi)-H_i**2)**(1/2) #equation 7
            V_i = (1/6)*pi*H_i*(3*A_i**2+H_i**2) #equation 6
        # this is actually wrong!
        #    R_i = (A_i**2-H_i**2)/(2*H_i) #equation 10
            R_i = S_i/(2*pi*H_i) #equation 10 re-do
        # this is wrong!
        #    theta_i = np.arcsin(A_i/R_i) #equation 12
            theta_i = np.arccos(1-S_i/(2*pi*R_i**2)) #equation 12 re-do
            V_b = V_d+V_i # equation 5
            R_d_list = []
            #for V in V_b:
            #    R_d_list.append(cube_root((3*V)/(pi*(2+np.cos(theta_d))*(1-np.cos(theta_d))**2)))
            R_d = ((3*V_b)/(pi*(2+np.cos(theta_d))*(1-np.cos(theta_d))**2))**(1/3) #equation 4
            # R_d = np.array(R_d_list)
            S_d = 2*pi*R_d**2*(1-np.cos(theta_d)) #equation 11
            A_d = R_d*np.sin(theta_d) #equation 3
            S_m = pi*(A_d**2-A_i**2) #equation 9
            S_o = pi*(A_t**2-A_d**2) #equation 2
            membrane_tension = sigma*(S_o+S_m+S_i)
            bending_energy = (K/2)*S_i*((2/R_i)-(2*C_i))**2
            surface_tension = gamma_d*S_d
            wetting_energy = gamma_m*S_m+gamma_i*S_i
            droplet_energy = surface_tension - wetting_energy
            turgor_pressure = P*V_i
            F = membrane_tension+bending_energy+surface_tension-wetting_energy+turgor_pressure #equation 1
            non_droplet_energy = membrane_tension+bending_energy+turgor_pressure
        
            outputs = {'theta_d':theta_d,'A_i':A_i,'V_i':V_i, 'R_i':R_i,
                       '2÷R_i':2/R_i,'theta_i':theta_i,'V_b':V_b,'R_d':R_d,
                      'S_d':S_d,'A_d':A_d,'S_m':S_m,'S_o':S_o,
                       'F':F,'membrane_tension':membrane_tension,
                       'bending_energy':bending_energy, 'bending_energy_flipped':-bending_energy,
                       'surface_tension':surface_tension,
                       'wetting_energy':-wetting_energy, 'wetting_energy_flipped':wetting_energy,
                       'droplet_energy':droplet_energy, 'turgor_pressure':turgor_pressure,
                       'non_droplet_energy':non_droplet_energy}
            return outputs
        
        plt.plot(H_i_range/100,llps_model()['F'])
        plt.xlabel('h')
        plt.ylabel('F')
        plt.tight_layout()
        plt.savefig(fig_dir+'llps_model_test.png')
        
      • post-analytical version
        def llps_model_recbrt(h=h_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                       C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i,
                       V_d=V_d, P=P):
            
            alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
            eta = (S_i**(3/2))/(6*pi**(1/2))
            C = sigma*pi*A_t**2+S_i*(gamma_m-gamma_i)
            
            F = ((sigma*alpha+(sigma-gamma_m)*beta+gamma_d*epsilon)* np.cbrt(((S_i**(3/2))/(pi**(1/2)))*(3*h/2-h**3)+3*V_d)**2+ (sigma-gamma_m)*S_i*h**2+ 2*K*(C_i*S_i**(1/2)-2*pi**(1/2)*h)**2+ P*eta*(3*h-2*h**3)+ C )
        
        #    print(alpha, beta, epsilon, eta, C)
        
            return F
        
        plt.plot(h_range,llps_model_recbrt())
        plt.xlabel('h')
        plt.ylabel('F')
        plt.tight_layout()
        plt.savefig(fig_dir+'llps_model_recbrt_test.png')
        
    • general plotting function
      • figure out some way to get the right units for each thing being plotted
      def plot_simulation(x=H_i_range/vesicle_diameter, output='F',
                       H_i=H_i_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                       C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P,
                       label=None, variable_df=variable_df):
          outputs = llps_model(H_i=H_i_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                    C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
          y = outputs[output]
          plt.plot(x,y,label=label)
          plt.ylabel(variable_df.loc[output]['description']+' ('+variable_df.loc[output]['units']+')')
      
      • testing
        plot_simulation()
        
    • normalized plots
      def plot_simulation_normalized(x=H_i_range/vesicle_diameter, output='F',
                      H_i=H_i_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                      C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P,
                                      label=None, color=None, linestyle='-', linewidth=None):
          outputs = llps_model(H_i=H_i_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                  C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d)
          y_sub = outputs[output]-np.min(outputs[output])
          y_norm = y_sub/np.max(y_sub)
          plt.plot(x,y_norm,label=label, color=color, linestyle=linestyle, linewidth=linewidth)
          plt.ylabel(output)
          plt.xlabel('h: invagination distance/vesicle diameter')
      
    • delta plots
      def plot_simulation_delta(x=H_i_range/vesicle_diameter, output='F',
                       H_i=H_i_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                       C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P,
                       label=None, variable_df=variable_df, color=None, linestyle='-', linewidth=None):
          outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                    C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
          output_0 = llps_model(H_i=H_i[0], sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                    C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
          y = outputs[output]-output_0[output]
          plt.plot(x,y,label=label, color=color,linestyle=linestyle, linewidth=linewidth)
          #plt.ylabel('∆ '+variable_df.loc[output]['description']+' ('+variable_df.loc[output]['units']+')')
          plt.xlabel('h: invagination distance/vesicle diameter')
      
    • show geometry for given invagination distance
      def show_geometry(H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                        C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P,
                        fig_width=10, fig_height=5,
                        num_plots=1, top=20, bottom=-250, left=-500, right=500,
                        intermediate=False):
          outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                   C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
      
          A_d = outputs['A_d']
          A_i = outputs['A_i']
          theta_d = outputs['theta_d']
          theta_i = outputs['theta_i']
          R_d = outputs['R_d']
          R_i = outputs['R_i']
      
          #fig, ax = plt.subplots(1,num_plots)
          #plt.subplot(1, num_plots, 1)
      
          fig = plt.figure()
          fig.set_size_inches((fig_width,fig_height))
          ax = fig.add_subplot(1, num_plots, 1)
      
          plt.hlines([0,0],[-A_t,A_d],[-A_d,A_t], color='r', label='$S_o$')
          plt.hlines([0,0],[-A_d,A_i],[-A_i,A_d], color='g', label='$S_m$')
      
          if V_d > 0:
            points_d = np.linspace(pi/2-theta_d, theta_d+pi/2, 1000)
            x_d = R_d*np.cos(points_d)
            y_d = R_d*np.sin(points_d)
            y_d_adj = y_d-np.min(y_d)
            plt.plot(x_d,-y_d_adj, color='b', label='$S_d$')
      
          points_i = np.linspace(pi/2-theta_i, theta_i+pi/2, 1000)
          x_i = R_i*np.cos(points_i)
          y_i = R_i*np.sin(points_i)
          y_i_adj = y_i-np.min(y_i)
          plt.plot(x_i,-y_i_adj, color='purple', label='$S_i$')
      
          ax.set_ylim(top=top, bottom=bottom)
          ax.set_xlim(left=left, right=right)
          ax.set_aspect('equal')
          ax.set(xlabel = 'distance along membrane (nm)',
                 ylabel = 'distance into cytoplasm (nm)')
          plt.legend()
          plt.tight_layout()
      
          if intermediate==True:
                 return fig, ax
      
      • testing
        show_geometry(H_i=100)
        
    • get geometry to feed into continuum model
      def get_geometry(H_i, N, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                        C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P):
          outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                   C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
      
          r = np.linspace(0, A_t, 3*N+1)
          z = np.full(3*N+1,0)
      
          A_d = outputs['A_d']
          A_i = outputs['A_i']
          theta_d = outputs['theta_d']
          theta_i = outputs['theta_i']
          R_d = outputs['R_d']
          R_i = outputs['R_i']
      
          r[N:2*N+1] = np.linspace(A_i,A_d,N+1)
          r[2*N:3*N+1] = np.linspace(A_d,A_t,N+1)
      
          points_i = np.linspace(pi/2, pi/2-theta_i, N+1)
          x_i = R_i*np.cos(points_i)
          y_i = R_i*np.sin(points_i)
          y_i_adj = y_i-np.min(y_i)
          r[0:N+1] = x_i
          z[0:N+1] = -y_i_adj
      
          return(r,z)
      
      • testing
        show_geometry(H_i=100)
        
    • 3D geometry
      • points
        %gui qt
        viewer = napari.Viewer(ndisplay=3)
        
        spherecoords = []
        for H_i in np.linspace(1,100,100):
            outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K, C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
        
            A_d = outputs['A_d']
            A_i = outputs['A_i']
            theta_d = outputs['theta_d']
            theta_i = outputs['theta_i']
            R_d = outputs['R_d']
            R_i = outputs['R_i']
            theta, phi = np.linspace(0,2* np.pi, 200), np.linspace(pi-theta_d, np.pi, 200)
            THETA, PHI = np.meshgrid(theta, phi)
            R = R_d
            X = R * np.sin(PHI) * np.cos(THETA)
            Y = R * np.sin(PHI) * np.sin(THETA)
            Z = R * np.cos(PHI)
            Z = Z-np.max(Z)
            height = np.full_like(Z, H_i)
            spherecoords.append(np.array([np.ravel(height), np.ravel(Z), np.ravel(Y), np.ravel(X)]).T)
        
        spherecoords = np.concatenate(spherecoords)
        viewer.add_points(np.array(spherecoords), name='droplet', size=20, face_color='cyan',
                          edge_color = 'transparent', blending='additive', opacity=0.6)
        
        spherecoords = []
        for H_i in np.linspace(1,100,100):
            outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K, C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
        
            A_d = outputs['A_d']
            A_i = outputs['A_i']
            theta_d = outputs['theta_d']
            theta_i = outputs['theta_i']
            R_d = outputs['R_d']
            R_i = outputs['R_i']
            theta, phi = np.linspace(0,2* np.pi, 200), np.linspace(pi-theta_i, np.pi, 200)
            THETA, PHI = np.meshgrid(theta, phi)
            R = R_i
            X = R * np.sin(PHI) * np.cos(THETA)
            Y = R * np.sin(PHI) * np.sin(THETA)
            Z = R * np.cos(PHI)
            Z = Z-np.max(Z)
            height = np.full_like(Z, H_i)
            spherecoords.append(np.array([np.ravel(height), np.ravel(Z), np.ravel(Y), np.ravel(X)]).T)
        
        spherecoords = np.concatenate(spherecoords)
        viewer.add_points(np.array(spherecoords), name='coat', size=20, face_color='purple',
                          edge_color = 'transparent', blending='additive', opacity=0.6)
        
        spherecoords = []
        for H_i in np.linspace(1,100,100):
            outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K, C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
        
            A_d = outputs['A_d']
            A_i = outputs['A_i']
            theta_d = outputs['theta_d']
            theta_i = outputs['theta_i']
            theta_m = 0.00000001
            R_d = outputs['R_d']
            R_i = outputs['R_i']
            theta, phi = np.linspace(0,2* np.pi, 200), np.linspace(pi-theta_m, np.pi, 200)
            THETA, PHI = np.meshgrid(theta, phi)
            R = 0.1*A_t/np.sin(theta_m)
            X = R * np.sin(PHI) * np.cos(THETA)
            Y = R * np.sin(PHI) * np.sin(THETA)
            Z = R * np.cos(PHI)
            Z = Z-np.max(Z)+0
            height = np.full_like(Z, H_i)
            spherecoords.append(np.array([np.ravel(height), np.ravel(Z), np.ravel(Y), np.ravel(X)]).T)
        
        spherecoords = np.concatenate(spherecoords)
        viewer.add_points(np.array(spherecoords), name='cell wall', size=20, face_color='red',
                          edge_color = 'transparent', blending='additive', opacity=0.6)
        
        
        • I don’t understand why this doesn’t really work
        vertices = spherecoords
        #faces = np.array([[0, 0, 1, 2], [0, 1, 2, 3], [0, 2, 3, 4]])
        faces_list = []
        counter = 0
        for H_i in range(len(H_i_range)):
            for point in range(38):
                pointa = counter + point
                faces_list.append([pointa, pointa+1, pointa+2])
            counter += 40
        faces = np.array(faces_list)
        values = np.linspace(0, 1, len(vertices))
        surface = (vertices, faces, values)
        
        viewer.add_surface(surface, colormap='viridis')  # add the surface
        
      • surface

        examples:

        vertices = np.array([[0, 0, 0, 0], [0, 10, 0, 20], [0, 20, 10, 0], [0, 0, 10, 10],
                             [1, 20, 0, 0], [1, 10, 0, 20], [1, 20, 10, 0], [1, 0, 10, 10]])
        faces = np.array([[0, 1, 2], [1, 2, 3],
                          [4,5,6],[5,6,7]])
        values = np.linspace(0.9, 1, len(vertices))
        surface = (vertices, faces, values)
        
        viewer = napari.view_surface(surface)  # add the surface
        
        facesloop = []
        
        for t in (0,1):
            vert = vertices[vertices[:,0] == t]
            facest = scipy.spatial.Delaunay(vert[:,2:]).simplices + len(vert)*t
            facesloop.append(facest)
        
        facesd = np.concatenate(facesloop)
        
        viewer = napari.view_surface((vertices, facesd, values))
        
        coords = np.array([np.ravel(height), np.ravel(Z), np.ravel(Y), np.ravel(X)]).T[:,1:]
        faces = scipy.spatial.Delaunay( coords[:,1:]).simplices
        values = np.linspace(0,1,len(coords))
        viewer = napari.view_surface((coords[:,1:], faces, values))
        
        %gui qt
        viewer = napari.Viewer(ndisplay=3)
        
        spherecoords = []
        faces = []
        for H_i in range(1,101):
            outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K, C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
        
            A_d = outputs['A_d']
            A_i = outputs['A_i']
            theta_d = outputs['theta_d']
            theta_i = outputs['theta_i']
            R_d = outputs['R_d']
            R_i = outputs['R_i']
            theta, phi = np.linspace(0,2* np.pi, 50), np.linspace(pi-theta_d, np.pi, 50)
            THETA, PHI = np.meshgrid(theta, phi)
            R = R_d
            X = R * np.sin(PHI) * np.cos(THETA)
            Y = R * np.sin(PHI) * np.sin(THETA)
            Z = R * np.cos(PHI)
            Z = Z-np.max(Z)
            height = np.full_like(Z, H_i)
            vertices = np.array([np.ravel(height), np.ravel(Z), np.ravel(Y), np.ravel(X)]).T
            spherecoords.append(vertices)
        
            #face = scipy.spatial.Delaunay(vertices[:,2:]).simplices + len(vertices)*(H_i-1)
            if H_i == 1:
                face = scipy.spatial.Delaunay(vertices[:,2:]).simplices + len(vertices)*(H_i-1)
            else:
                face = faces[-1] + len(vertices)
            faces.append(face)
        
        spherecoords = np.concatenate(spherecoords)
        faces = np.concatenate(faces)
        surface_values = np.full(len(spherecoords),1)
        viewer.add_surface((spherecoords,faces,surface_values), name='droplet', colormap='green',
                           blending='translucent', shading='smooth')
        
        spherecoords = []
        faces = []
        for H_i in range(1,101):
            outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K, C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
        
            A_d = outputs['A_d']
            A_i = outputs['A_i']
            theta_d = outputs['theta_d']
            theta_i = outputs['theta_i']
            R_d = outputs['R_d']
            R_i = outputs['R_i']
            theta, phi = np.linspace(0,2* np.pi, 50), np.linspace(pi-theta_i, np.pi, 50)
            THETA, PHI = np.meshgrid(theta, phi)
            R = R_i
            X = R * np.sin(PHI) * np.cos(THETA)
            Y = R * np.sin(PHI) * np.sin(THETA)
            Z = R * np.cos(PHI)
            Z = Z-np.max(Z)
            height = np.full_like(Z, H_i)
            vertices = np.array([np.ravel(height), np.ravel(Z), np.ravel(Y), np.ravel(X)]).T
            spherecoords.append(vertices)
        
            if H_i == 1:
                face = scipy.spatial.Delaunay(vertices[:,2:]).simplices + len(vertices)*(H_i-1)
            else:
                face = faces[-1] + len(vertices)
            faces.append(face)
        
        spherecoords = np.concatenate(spherecoords)
        faces = np.concatenate(faces)
        surface_values = np.full(len(spherecoords),1)
        viewer.add_surface((spherecoords,faces,surface_values), name='coat', colormap='bop blue',
                           blending='additive', shading='smooth')
        
        spherecoords = []
        faces = []
        for H_i in range(1,101):
            outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K, C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
        
            A_d = outputs['A_d']
            A_i = outputs['A_i']
            theta_d = outputs['theta_d']
            theta_i = outputs['theta_i']
            theta_m = 0.00000001
            R_d = outputs['R_d']
            R_i = outputs['R_i']
            theta, phi = np.linspace(0,2* np.pi, 50), np.linspace(pi-theta_m, np.pi, 50)
            THETA, PHI = np.meshgrid(theta, phi)
            R = 0.1*A_t/np.sin(theta_m)
            X = R * np.sin(PHI) * np.cos(THETA)
            Y = R * np.sin(PHI) * np.sin(THETA)
            Z = R * np.cos(PHI)
            Z = Z-np.max(Z)+0
            height = np.full_like(Z, H_i)
            vertices = np.array([np.ravel(height), np.ravel(Z), np.ravel(Y), np.ravel(X)]).T
            spherecoords.append(vertices)
        
            #face = scipy.spatial.Delaunay(vertices[:,2:]).simplices + len(vertices)*(H_i-1)
            if H_i == 1:
                face = scipy.spatial.Delaunay(vertices[:,2:]).simplices + len(vertices)*(H_i-1)
            else:
                face = faces[-1] + len(vertices)
            faces.append(face)
        
        spherecoords = np.concatenate(spherecoords)
        faces = np.concatenate(faces)
        surface_values = np.full(len(spherecoords),1)
        viewer.add_surface((spherecoords,faces,surface_values), name='cell wall', colormap='bop blue',
                           blending='opaque', shading='smooth')
        
        def update_frame():
            """Update frame or time"""
            viewer.text_overlay.text = 'h = %1.2f\nblue: membrane\ngreen: condensate' % ((viewer.dims.current_step[0]+1)/100)
        
        viewer.dims.events.current_step.connect(lambda event: update_frame())
        
        viewer.text_overlay.visible=True
        viewer.text_overlay.font_size=20
        viewer.text_overlay.color = 'white'
        viewer.scale_bar.unit = 'nm'
        viewer.scale_bar.visible=True
        viewer.scale_bar.font_size=10
        

        animation script

        from napari_animation import Animation
        
        animation = Animation(viewer)
        viewer.reset_view()
        viewer.dims.current_step = (0, 0,0,0)
        animation.capture_keyframe(steps=15)
        animation.capture_keyframe(steps=15)
        viewer.camera.zoom = 2.6
        viewer.camera.angles=(-0.026467524771271397, 0.5719902891444765, 11.550928988495166)
        animation.capture_keyframe(steps=50)
        viewer.dims.current_step = (viewer.dims.range[0][1]-1, 0,0,0)
        animation.capture_keyframe(steps=int(viewer.dims.range[0][1]-1))
        viewer.reset_view()
        animation.capture_keyframe(steps=50)
        viewer.dims.current_step = (0, 0,0,0)
        animation.capture_keyframe(steps=int(viewer.dims.range[0][1]-1))
        animation.animate('~/google_drive/grad_school/db_lab/code/analysis/llps_endocytosis_modeling/scripted_surface_animation.mov')
        
      • probably don’t need this stuff
        num_plots = 1
        outputs = llps_model(H_i=90, sigma=sigma, A_t=A_t, S_i=S_i, K=K, C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
        
        A_d = outputs['A_d']
        A_i = outputs['A_i']
        theta_d = outputs['theta_d']
        theta_i = outputs['theta_i']
        R_d = outputs['R_d']
        R_i = outputs['R_i']
        
        #fig, ax = plt.subplots(1,num_plots)
        #plt.subplot(1, num_plots, 1)
        
        #fig = plt.figure()
        #fig.set_size_inches((fig_width,fig_height))
        #ax = fig.add_subplot(1, num_plots, 1)
        
        #plt.hlines([0,0],[-A_t,A_d],[-A_d,A_t], color='r', label='$S_o$')
        #plt.hlines([0,0],[-A_d,A_i],[-A_i,A_d], color='g', label='$S_m$')
        
        if V_d > 0:
            points_d = np.linspace(pi/2-theta_d, theta_d+pi/2, 1000)
            x_d = R_d*np.cos(points_d)
            y_d = R_d*np.sin(points_d)
            y_d_adj = y_d-np.min(y_d)
        
        points_i = np.linspace(pi/2-theta_i, theta_i+pi/2, 1000)
        x_i = R_i*np.cos(points_i)
        y_i = R_i*np.sin(points_i)
        y_i_adj = y_i-np.min(y_i)
        
        xy_i = np.array([y_i,x_i]).T
        viewer = napari.Viewer()
        viewer.add_points(xy_i)
        
        nx, ny, nz = (3, 3, 4)
        x = np.linspace(0, 10, nx)
        y = np.linspace(0, 10, ny)
        z = np.linspace(0, 10, nz)
        
        xv, yv, zv = np.meshgrid(x, y, z)
        zyxv = np.array([np.ravel(zv), np.ravel(yv), np.ravel(xv)]).T
        zyxv
        
        viewer = napari.Viewer()
        viewer.add_points(zyxv)
        
        import numpy as np
        import matplotlib.pyplot as plt
        from mpl_toolkits import mplot3d
        
        ax = plt.axes(projection="3d")
        xlist=np.linspace(-1.0,1.0,10)*100
        ylist=np.linspace(-1.0,1.0,10)*100
        r=np.linspace(1.0,1.0,10)*100
        X,Y= np.meshgrid(xlist,ylist)
        
        Z=np.sqrt(r**2-X**2-Y**2) #Use np.sqrt like you had before
        
        spherecoords = np.array([np.ravel(Z), np.ravel(Y), np.ravel(X)]).T
        viewer.add_points(spherecoords)
        
        vertices = spherecoords
        faces = np.array([[0, 50, 99], [1, 2, 3]])
        values = np.linspace(1, 1, len(vertices))
        surface = (vertices, faces, values)
        
        viewer.add_surface(surface)
        
        cp=ax.plot_wireframe(X,Y,Z,color="r")
        cp=ax.plot_wireframe(X,Y,-Z,color="r") # Now plot the bottom half
        
        
        plt.title('2D Contour Plot of The unit sphere')
        plt.show()
        
    • side-by-side delta plot and equilibrium geometry
      def plot_simulation_geometry(x=H_i_range/vesicle_diameter, output='F',
                       H_i=H_i_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                       C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P,
                       fig_width=15, fig_height=5, label=None, variable_df=variable_df,
                      num_plots=2, top=20, bottom=-250, left=-500, right=500):
          outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                    C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
          min_energy_ix = np.argmin(outputs['F'])
          H_i_eq = H_i[min_energy_ix]
          fig, ax1 = show_geometry(H_i=H_i_eq, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                         C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P,
                         fig_width=fig_width, fig_height=fig_height,
                         num_plots=num_plots, top=top, bottom=bottom, left=left, right=right,
                         intermediate=True)
          plt.title('geometry at equilibrium')
          ax2 = fig.add_subplot(1,2,2)
          #y = outputs[output]
          #plt.plot(x,y,label=label)
          #plt.ylabel(output)
      
          plot_simulation_delta(x=x, output=output,
                       H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                       C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P,
                       label=label, variable_df=variable_df)
      
          plt.tight_layout()
      
      • testing
        plot_simulation_geometry()
        
    • animated side-by-side geometry and energy plot, scan volume
      V_d_range=np.linspace(V_d*0.1,V_d*10,100)
      def animate_simulation_geometry(x=H_i_range/vesicle_diameter, output='F',
                      H_i=H_i_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                      C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d_range=V_d_range, P=P,
                      fig_width=15, fig_height=5, label=None, color=None, variable_df=variable_df,
                      geometry_top=20, geometry_bottom=-250, geometry_left=-500, geometry_right=500,
                      plot_top=20, plot_bottom=-250, plot_left=-500, plot_right=500):
          i = 0
          V_d=V_d_range[i]
      
          fig, (ax_geometry, ax_plot) = plt.subplots(nrows=1, ncols=2, figsize=(15,5))
      
          outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                      C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
          output_0 = llps_model(H_i=H_i[0], sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                  C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
          min_energy_ix = np.argmin(outputs['F'])
          H_i_eq = H_i[min_energy_ix]
          outputs_eq = llps_model(H_i=H_i_eq, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                  C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
      
          A_d = outputs_eq['A_d']
          A_i = outputs_eq['A_i']
          theta_d = outputs_eq['theta_d']
          theta_i = outputs_eq['theta_i']
          R_d = outputs_eq['R_d']
          R_i = outputs_eq['R_i']
      
          top = geometry_top
          bottom = geometry_bottom
          left = geometry_left
          right = geometry_right
      
          S_o_line_l, = ax_geometry.plot([-A_t,-A_d],[0,0], color='b', label='$S_o$')
          S_o_line_r, = ax_geometry.plot([A_d,A_t],[0,0], color='b')
          S_m_line_l, = ax_geometry.plot([-A_d,-A_i],[0,0], color='cyan', label='$S_m$')
          S_m_line_r, = ax_geometry.plot([A_i,A_d],[0,0], color='cyan')
          if V_d > 0:
              points_d = np.linspace(pi/2-theta_d, theta_d+pi/2, 1000)
              x_d = R_d*np.cos(points_d)
              y_d = R_d*np.sin(points_d)
              y_d_adj = y_d-np.min(y_d)
              S_d_line, = ax_geometry.plot(x_d,-y_d_adj, color='g', label='$S_d$')
              if geometry_bottom == None:
                  bottom = y_d_adj.max()*(-1.1)
              if geometry_left == None:
                  left = x_d.max()*(-1.2)
              if geometry_right == None:
                  right = x_d.max()*(1.2)
          points_i = np.linspace(pi/2-theta_i, theta_i+pi/2, 1000)
          x_i = R_i*np.cos(points_i)
          y_i = R_i*np.sin(points_i)
          y_i_adj = y_i-np.min(y_i)
          S_i_line, = ax_geometry.plot(x_i,-y_i_adj, color='purple', label='$S_i$')
      
          ax_geometry.set_ylim(top=top, bottom=bottom)
          ax_geometry.set_xlim(left=left, right=right)
          ax_geometry.set_aspect('equal')
          ax_geometry.set(xlabel = 'distance along membrane (nm)',
                  ylabel = 'distance into cytoplasm (nm)')
          ax_geometry.legend()
          ax_geometry.set_title('geometry at equilibrium')
      
      
          surface_tension = outputs['surface_tension']-output_0['surface_tension']
          wetting_energy = outputs['wetting_energy']-output_0['wetting_energy']
          fd = outputs['droplet_energy']-output_0['droplet_energy']
          plot_st, = ax_plot.plot(x,surface_tension,label='$\gamma_d S_d$', color='green')
          plot_we, = ax_plot.plot(x,wetting_energy,label='$\gamma_m S_m$', color='cyan')
          plot_fd, = ax_plot.plot(x,fd,label='$F_d$', color='black')
          plt.hlines(0,0,1,linestyles=':',colors='gray')
          if plot_top == None:
              ax_plot.set_ylim(top=surface_tension.max()*1.2,
                              bottom=wetting_energy.min()-surface_tension.max()*0.2)
          else:
              ax_plot.set_ylim(top=plot_top, bottom=plot_bottom)
          ax_plot.legend()
          ax_plot.set(xlabel='h: invagination distance/vesicle diameter',
              ylabel='∆ energy $(pN \cdot nm)$')
      
          plt.tight_layout()
      
          def animate(i):
              V_d=V_d_range[i]
      
              outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                          C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
              output_0 = llps_model(H_i=H_i[0], sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                      C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
              min_energy_ix = np.argmin(outputs['F'])
              H_i_eq = H_i[min_energy_ix]
              outputs_eq = llps_model(H_i=H_i_eq, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                      C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
      
              A_d = outputs_eq['A_d']
              A_i = outputs_eq['A_i']
              theta_d = outputs_eq['theta_d']
              theta_i = outputs_eq['theta_i']
              R_d = outputs_eq['R_d']
              R_i = outputs_eq['R_i']
      
              top = geometry_top
              bottom = geometry_bottom
              left = geometry_left
              right = geometry_right
      
              #S_o_line = ax_geometry.hlines([0,0],[-A_t,A_d],[-A_d,A_t], color='r', label='$S_o$')
              S_o_line_l.set_data([-A_t,-A_d],[0,0])
              S_o_line_r.set_data([A_d,A_t],[0,0])
              #S_m_line = ax_geometry.hlines([0,0],[-A_d,A_i],[-A_i,A_d], color='g', label='$S_m$')
              S_m_line_l.set_data([-A_d,-A_i],[0,0])
              S_m_line_r.set_data([A_i,A_d],[0,0])
              if V_d > 0:
                  points_d = np.linspace(pi/2-theta_d, theta_d+pi/2, 1000)
                  x_d = R_d*np.cos(points_d)
                  y_d = R_d*np.sin(points_d)
                  y_d_adj = y_d-np.min(y_d)
                  S_d_line.set_data(x_d,-y_d_adj)
                  if geometry_bottom == None:
                      bottom = y_d_adj.max()*(-1.1)
                  if geometry_left == None:
                      left = x_d.max()*(-1.2)
                  if geometry_bottom == None:
                      right = x_d.max()*(1.2)
              points_i = np.linspace(pi/2-theta_i, theta_i+pi/2, 1000)
              x_i = R_i*np.cos(points_i)
              y_i = R_i*np.sin(points_i)
              y_i_adj = y_i-np.min(y_i)
              S_i_line.set_data(x_i,-y_i_adj)
              ax_geometry.set_ylim(top=top, bottom=bottom)
              ax_geometry.set_xlim(left=left, right=right)
              ax_geometry.set_aspect('equal')
      
              surface_tension = outputs['surface_tension']-output_0['surface_tension']
              wetting_energy = outputs['wetting_energy']-output_0['wetting_energy']
              fd = outputs['droplet_energy']-output_0['droplet_energy']
              plot_st.set_data(x,surface_tension)
              plot_we.set_data(x,wetting_energy)
              plot_fd.set_data(x,fd)
              if plot_top == None:
                  ax_plot.set_ylim(top=surface_tension.max()*1.2,
                                  bottom=wetting_energy.min()-surface_tension.max()*0.2)
      
              plt.tight_layout()
          ani = animation.FuncAnimation(
              fig, animate, frames=len(V_d_range), interval=100)
          ani.save(os.path.join(animation_dir, 'vdscan.mp4'), dpi=300)
      
          plt.show()
      
    • side-by-side F plots and parameter phase diagram
      def boundary_test(x_param, y_param, x_range, boundary, params=params,
                        factors=[0.7,1,1.3], y_offset=0.2, x_scale='log'):
          fparams = params.copy()
          fig, (ax1, ax2)  = plt.subplots(nrows=1, ncols=2, figsize=(10,4))
      
          if x_scale == 'log':
              ax1.semilogx(x_range, boundary, color='black')
          else:
              ax1.plot(x_range, boundary, color='black')
      
          for factor in factors:
              middle = len(x_range)//2
              x_pick = x_range[middle]*factor
              boundary_pick = boundary[middle]+(factor-1)*np.sqrt(boundary[middle]**2)*y_offset
              fparams[x_param] = x_pick
              fparams[y_param] = boundary_pick
              F = llps_model_recbrt(h=fparams['h_range'], sigma=fparams['sigma'], A_t=fparams['A_t'],
                      S_i=fparams['S_i'], K=fparams['K'], C_i=fparams['C_i'], gamma_d=fparams['gamma_d'],
                      gamma_m=fparams['gamma_m'], gamma_i=fparams['gamma_i'], V_d=fparams['V_d'], P=fparams['P'])
              delta_F = F-F[0]
      
              ax1.scatter(x_pick, boundary_pick)
              ax2.plot(h_range,delta_F)
      
          ax2.hlines(0, 0,1.05, color='black')
      
          ax1.set(xlabel=x_param, ylabel=y_param)
          ax2.set(xlabel='h', ylabel='∆F')
          plt.tight_layout()
      
    • global stability parameter analysis [4/4]
      • DONE \(K\)

        from therm_spont_k

        def thermo_spont_K(sigma=sigma, gamma_m=gamma_m, gamma_d=gamma_d,
                          V_d=V_d, S_i=S_i, P=P, C_i=C_i):
            
            alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
            eta = (S_i**(3/2))/(6*pi**(1/2))
            
            K = ((-(sigma*alpha+(sigma-gamma_m)*beta+gamma_d*epsilon)*
                  (np.cbrt(3*V_d+(S_i**(3/2))/(2*pi**(1/2)))**2-
                    np.cbrt(3*V_d)**2)+2*P*eta-(sigma-gamma_m)*S_i-
                  3*P*eta)/(8*pi-8*C_i*(pi*S_i)**(1/2)))
            
            return K
        
        • testing
          V_d_range = np.logspace(6,9)
          K_boundary = thermo_spont_K(V_d=V_d_range)
          boundary_test(x_param='V_d', y_param='K', x_range=V_d_range, boundary=K_boundary)
          
          C_i_range = np.logspace(-4,-2.1,1000)
          K_boundary = thermo_spont_K(C_i=C_i_range)
          boundary_test(x_param='C_i', y_param='K', x_range=C_i_range, boundary=K_boundary,
                        y_offset=-0.3)
          
          • this kind of makes sense - the lower the spontaneous curvature (flatter), the more the energetics are sensitive to rigidity
          x_range = np.logspace(6,9)
          thermo_boundary = thermo_spont_K(V_d=x_range)
          kinetic_boundary = kinetic_spont_K(V_d=x_range)
          both_boundary_test(x_param='V_d', y_param='K', x_range=x_range,
                             thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                             y_offset=2)
          
      • DONE \(P\)
        def thermo_spont_P(sigma=sigma, gamma_m=gamma_m, gamma_d=gamma_d,
                          V_d=V_d, S_i=S_i, K=K, C_i=C_i):
            
            alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
            eta = (S_i**(3/2))/(6*pi**(1/2))
            
            P = ((sigma*alpha + (sigma - gamma_m)*beta + gamma_d*epsilon)*(np.cbrt(3*V_d)**2 -
                 np.cbrt(S_i**(3/2)/(2*np.sqrt(pi)) + 3*V_d)**2) -
                 (sigma-gamma_m)*S_i - 8*K*(pi - C_i*np.sqrt(pi*S_i)))/ eta
            
            return P
        
        • testing
          x_range = np.logspace(6,9)
          boundary = thermo_spont_P(V_d=x_range)
          boundary_test(x_param='V_d', y_param='P', x_range=x_range, boundary=boundary,
                        y_offset=-0.3)
          
      • DONE \(\sigma\)
        def thermo_spont_sigma(P=P, gamma_m=gamma_m, gamma_d=gamma_d,
                          V_d=V_d, S_i=S_i, K=K, C_i=C_i):
            
            alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
            eta = (S_i**(3/2))/(6*pi**(1/2))
        
            sigma = ((-gamma_m*S_i+8*K*(pi-C_i*np.sqrt(pi*S_i))+P*eta)/
                    ((3*V_d)**(2/3)-((S_i**(3/2))/
                        (2*np.sqrt(pi))+3*V_d)**(2/3))-
                     gamma_d*epsilon+gamma_m*beta)/(alpha+beta-
                    (S_i)/((3*V_d)**(2/3)-((S_i**(3/2))/
                            (2*np.sqrt(pi))+3*V_d)**(2/3)))
            
            return sigma
        
        • testing
          x_range = np.logspace(6,9)
          boundary = thermo_spont_sigma(V_d=x_range)
          boundary_test(x_param='V_d', y_param='sigma', x_range=x_range, boundary=boundary,
                        y_offset=-0.02)
          
      • DONE \(C_{i}\)
        def thermo_spont_C_i(P=P, gamma_m=gamma_m, gamma_d=gamma_d,
                          V_d=V_d, S_i=S_i, K=K, sigma=sigma):
            
            alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
            eta = (S_i**(3/2))/(6*pi**(1/2))
        
            C_i = (pi-((sigma*alpha+(sigma-gamma_m)*beta+gamma_d*epsilon)*(np.cbrt(3*V_d)**2-np.cbrt((S_i**(3/2))/(2*np.sqrt(pi))+3*V_d)**2)-(sigma-gamma_m)*S_i-P*eta)/(8*K))/(np.sqrt(pi*S_i))
            
            return C_i
        
        • testing
          x_range = np.logspace(2,7)
          boundary = thermo_spont_C_i(V_d=x_range)
          boundary_test(x_param='V_d', y_param='C_i', x_range=x_range, boundary=boundary,
                        y_offset=0.2)
          
          • it doesn’t make sense that C_i would ever be negative…
    • local stability parameter analysis [8/8]
      • DONE \(K\)

        from kinetic_spont_k

        def kinetic_spont_K(sigma=sigma, gamma_m=gamma_m, gamma_d=gamma_d,
                          V_d=V_d, S_i=S_i, P=P, C_i=C_i):
            
            alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
            eta = (S_i**(3/2))/(6*pi**(1/2))
            
            K = (((sigma*alpha+(sigma-gamma_m)*beta+gamma_d*epsilon)*S_i**(3/2))/(np.sqrt(pi)*np.cbrt(3*V_d)) + 3*P*eta)/(8*C_i*np.sqrt(pi*S_i))
            
            return K
        
        • testing
          x_range = np.logspace(6,9)
          boundary = kinetic_spont_K(V_d=x_range)
          boundary_test(x_param='V_d', y_param='K', x_range=x_range, boundary=boundary,
                        y_offset=0.5)
          
      • DONE \(P\)

        from kinetic_spont_p

        def kinetic_spont_P(sigma=sigma, gamma_m=gamma_m, gamma_d=gamma_d,
                          V_d=V_d, S_i=S_i, K=K, C_i=C_i):
            
            alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
            eta = (S_i**(3/2))/(6*pi**(1/2))
            
            P = (-((sigma*alpha + (sigma - gamma_m)*beta +gamma_d*epsilon)*S_i**(3/2))/
                 (np.sqrt(pi)*np.cbrt(3*V_d)) + 8*K*C_i*np.sqrt(pi*S_i))/(3*eta)
            
            return P
        
        • testing
          x_range = np.logspace(6,9)
          boundary = kinetic_spont_P(V_d=x_range)
          boundary_test(x_param='V_d', y_param='P', x_range=x_range, boundary=boundary,
                        y_offset=-0.5)
          
      • DONE \(\sigma\)
        ---------------------------------------------------------------------------
        NameError                                 Traceback (most recent call last)
        <ipython-input-35-594c148ecf9a> in <module>
        ----> 1 def kinetic_spont_sigma(P=P, gamma_m=gamma_m, gamma_d=gamma_d,
              2                   V_d=V_d, S_i=S_i, K=K, C_i=C_i):
              3 
              4     alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
              5     beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
        
        NameError: name 'P' is not defined
        
        • testing
          x_range = np.logspace(6,9)
          boundary = kinetic_spont_sigma(V_d=x_range)
          boundary_test(x_param='V_d', y_param='sigma', x_range=x_range, boundary=boundary,
                        y_offset=-0.5)
          

          :END:

      • DONE \(C_{i}\)

        from kinetic_spont_ci

        def kinetic_spont_C_i(sigma=sigma, gamma_m=gamma_m, gamma_d=gamma_d,
                          V_d=V_d, S_i=S_i, K=K, P=P):
            
            alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
            eta = (S_i**(3/2))/(6*pi**(1/2))
            
            C_i = (((sigma*alpha + (sigma - gamma_m)*beta +gamma_d*epsilon)*S_i**(3/2))/
                 (np.sqrt(pi)*np.cbrt(3*V_d)) + 3*P*eta)/(8*K*np.sqrt(pi*S_i))
            
            return C_i
        
        • testing
          x_range = np.logspace(6,9)
          boundary = kinetic_spont_C_i(V_d=x_range)
          boundary_test(x_param='V_d', y_param='C_i', x_range=x_range, boundary=boundary,
                        y_offset=0.5)
          
      • CANCELLED \(\gamma_m\)
      • CANCELLED \(\gamma_d\)
      • CANCELLED \(S_{i}\)
      • DONE \(V_{d}\)

        from kinetic_spont_ci

        def kinetic_spont_V_d(sigma=sigma, gamma_m=gamma_m, gamma_d=gamma_d,
                          C_i=C_i, S_i=S_i, K=K, P=P):
            
            alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
            eta = (S_i**(3/2))/(6*pi**(1/2))
            
            V_d = (1/3)*(((sigma*alpha + (sigma - gamma_m)*beta +gamma_d*epsilon)*S_i**(3/2))/
                 ((8*K*C_i*np.sqrt(pi*S_i)-3*P*eta)*np.sqrt(pi)))**3
            
            return V_d
        
        • testing
          x_range = np.linspace(0.0002, 0.002)
          boundary = kinetic_spont_V_d(P=x_range)
          boundary_test(x_param='P', y_param='V_d', x_range=x_range, boundary=boundary,
                        y_offset=-5, x_scale='linear')
          
    • F’(0)
      def Fprime(h=0, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                 C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i,
                 V_d=V_d, P=P):
      
          alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
          beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
          epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
          eta = (S_i**(3/2))/(6*pi**(1/2))
          C = sigma*pi*A_t**2+S_i*(gamma_m-gamma_i)
      
          Fp = ( 3 * P * eta + ( -6 * ( h )**( 2 ) * P * eta + ( -8 * K * ( np.pi )**( 1/2) *
          C_i * ( S_i )**( 1/2 ) + ( ( 1 + -2 * ( h )**( 2 ) ) * (np.pi )**( -1/2 ) * (
          S_i)**( 3/2 ) * ( ( ( 3/2 * h + -1 * ( h )**( 3 ) ) * ( np.pi )**( -1/2 ) * (S_i
          )**( 3/2 ) + 3 * V_d ) )**( -1/3 ) * (sigma*alpha + ( epsilon * gamma_d + beta *
          ( sigma + -1 * gamma_m ) ) ) + 2 * h * (8 * K * np.pi + S_i * (sigma + -1 *
          gamma_m ) ) ) )) )
      
          return Fp
      
    • side-by-side F plots with both phase boundaries
      def both_boundary_test(x_param, y_param, x_range, thermo_boundary, kinetic_boundary,
                             params=params, factors=[0.7,1,1.3], y_offset=0.2, x_scale='log', y_scale='linear',
                             variable_df=variable_df, fill=False):
          fparams = params.copy()
          fig, (ax1, ax2)  = plt.subplots(nrows=1, ncols=2, figsize=(12,5))
          
      
          # if x_scale == 'log':
          #     ax1.semilogx(x_range, thermo_boundary, color='black', label='thermodynamic')
          #     ax1.semilogx(x_range, kinetic_boundary, color='blue', label='kinetic')
          # # else:
          #     ax1.plot(x_range, boundary, color='black')
          
          if y_offset > 0:
              fill_range = np.max(kinetic_boundary)
          if y_offset < 0:
              fill_range = np.min(kinetic_boundary)
              
          ax1.loglog(x_range, thermo_boundary, color='orange', label='global stability boundary')
          ax1.loglog(x_range, kinetic_boundary, color='blue', label='local stability boundary')
          if fill == True:
            ax1.fill_between(x_range, y1=thermo_boundary, y2=fill_range*1.3, color='yellow', alpha=0.8, label='globally stable')
            ax1.fill_between(x_range, y1=kinetic_boundary, y2=fill_range*1.3, color='blue', alpha=0.8, label='locally stable')
      
          ax1.legend()
      
          for i, factor in enumerate(factors):
              middle = len(x_range)//2
              x_pick = x_range[middle]*factor
              for j, boundary in enumerate([thermo_boundary, kinetic_boundary]):
                boundary_pick = boundary[middle]+(factor-1)*np.sqrt(thermo_boundary[middle]**2)*y_offset
                fparams[x_param] = x_pick
                fparams[y_param] = boundary_pick
                F = llps_model_recbrt(h=fparams['h_range'], sigma=fparams['sigma'], A_t=fparams['A_t'],
                        S_i=fparams['S_i'], K=fparams['K'], C_i=fparams['C_i'], gamma_d=fparams['gamma_d'],
                        gamma_m=fparams['gamma_m'], gamma_i=fparams['gamma_i'], V_d=fparams['V_d'], P=fparams['P'])
                delta_F = F-F[0]
      
                # ax1.scatter(x_pick, boundary_pick, linewidths=5, color=((j)/2, (j)/2, (j+1)/2, (i+1)/len(factors)))
                # ax2.plot(h_range,delta_F, color=((j)/2, (j)/2, (j+1)/2, (i+1)/len(factors)))
                
                ax1.scatter(x_pick, boundary_pick, linewidths=5)
                ax2.plot(h_range,delta_F)
      
          ax2.hlines(0, -0.05,1.05, color='black')
      
          ax1.set(xlabel=variable_df.loc[x_param]['description']+' ('+variable_df.loc[x_param]['units']+')',
                  ylabel=variable_df.loc[y_param]['description']+' ('+variable_df.loc[y_param]['units']+')')
      #            xscale=x_scale)
          ax1.set_xscale(x_scale)
          ax1.set_yscale(y_scale)
          ax2.set(xlabel='h: relative invagination height', ylabel='∆F ($pN \cdot nm$)')
          plt.tight_layout()
      

result plots [8/9]

  • plot all outputs vs. invagination height for default parameters
    • as individual plots
        outputs = llps_model()
        outputs.pop('theta_d')
        # outputs.pop('2÷R_i')
        variables = list(outputs.keys())
      
        for variable in variables:
            plt.figure()
            plot_simulation(output=variable)
            plt.xlabel('h: invagination distance/vesicle diameter')
            plt.savefig(fig_dir+'%s-vs-H_i.png' %(variable))
      
    • as subplots of one figure
        outputs = llps_model()
        outputs.pop('theta_d')
        variables = list(outputs.keys())
        nplots = len(variables)
        width = 3
        height = nplots//3+1
        plt.figure(figsize=(width*5,height*4))
      
        for i, variable in enumerate(variables):
            plt.subplot(height, width, i+1)
            plot_simulation(output=variable)
            plt.xlabel('h: invagination distance/vesicle diameter')
      
        plt.tight_layout()
        plt.savefig(fig_dir+'all_outputs.png')
      
  • change in all outputs vs. invagination height for default parameters
    • as individual plots
          outputs = llps_model()
          outputs.pop('theta_d')
          variables = outputs.keys()
      
          for variable in variables:
              plt.figure()
              plot_simulation_delta(output=variable)
              plt.tight_layout()
              plt.savefig(fig_dir+'%s-vs-H_i_delta.png' %(variable))
      
    • as subplots of one figure
        outputs = llps_model()
        outputs.pop('theta_d')
        variables = outputs.keys()
        nplots = len(variables)
        width = 3
        height = nplots//3+1
        plt.figure(figsize=(width*5,height*4))
      
        for i, variable in enumerate(variables):
            plt.subplot(height, width, i+1)
            plot_simulation_delta(output=variable)
      
        plt.tight_layout()
        plt.savefig(fig_dir+'all_outputs_delta.png')
      
  • normalized geometry terms
    plot_params = [
        ['2÷R_i',r'$\frac{1}{R_i}$ $(nm^{-1})$','purple','-'],
        ['V_i','$V_i$ $(nm^3)$','brown','-'],
        ['S_d','$S_d$ $(nm^2)$','green',':'],
        ['S_m','$S_m$ $(nm^2)$','cyan','-'],
        ['S_o','$S_o$ $(nm^2)$','blue','-']
    ]
    
    plt.figure(figsize=(6,6))
    
    for variable, label, color, linestyle in plot_params:
        plot_simulation_normalized(output=variable, H_i=H_i_range,
                        sigma=model_params['sigma'], A_t=model_params['A_t'],
                        S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                        gamma_d=model_params['gamma_d'],
                        gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                        gamma_i=model_params['gamma_1'], V_d=model_params['Vd'],
                        P=model_params['P'], label=label, color=color, linestyle=linestyle,
                        linewidth=4)
    
    plt.legend()
    #plt.ylim(-100, 100)
    plt.ylabel('normalized shape change')
    #plt.title(title)
    plt.tight_layout()
    plt.savefig(fig_dir+'normalized_shapes.png')
    plt.savefig(fig_dir+'normalized_shapes.svg')
    
  • change in surface areas
    plot_params = [
        ['S_d','$S_d$','green','-'],
        ['S_m','$S_m$','cyan','-']
    ]
    
    plt.figure(figsize=(6,4))
    
    for variable, label, color, linestyle in plot_params:
        plot_simulation_delta(output=variable,
                        sigma=model_params['sigma'], A_t=model_params['A_t'],
                        S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                        gamma_d=model_params['gamma_d'],
                        gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                        gamma_i=model_params['gamma_1'], V_d=model_params['Vd'],
                        P=model_params['P'], label=label, color=color, linestyle=linestyle,
                        linewidth=4)
    plt.ylabel('∆ surface area $(nm^2)$')
    plt.legend()
    plt.tight_layout()
    plt.savefig(fig_dir+'sd_si_vs-H_i_delta.svg')
    plt.savefig(fig_dir+'sd_si_vs-H_i_delta.png')
    
    plt.figure(figsize=(6,4))
    
    for variable, label, color, linestyle in plot_params:
        plot_simulation_delta(x=0.3*H_i_range/vesicle_diameter, H_i=0.3*H_i_range,output=variable,
                        sigma=model_params['sigma'], A_t=model_params['A_t'],
                        S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                        gamma_d=model_params['gamma_d'],
                        gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                        gamma_i=model_params['gamma_1'], V_d=model_params['Vd'],
                        P=model_params['P'], label=label, color=color, linestyle=linestyle,
                        linewidth=4)
    plt.ylabel('∆ surface area $(nm^2)$')
    plt.legend()
    plt.tight_layout()
    plt.savefig(fig_dir+'sd_si_vs-H_i_delta_zoom.svg')
    plt.savefig(fig_dir+'sd_si_vs-H_i_delta_zoom.png')
    
    plot_params = [
        ['surface_tension','$\gamma_d S_d$','green','-'],
        ['wetting_energy','$- \gamma_m S_m$','cyan','-'],
        ['droplet_energy','$F_d$','black','-']
    ]
    
    plt.figure(figsize=(6,4))
    plt.hlines(0,0,1,linestyles=':',colors='gray')
    
    for variable, label, color, linestyle in plot_params:
        plot_simulation_delta(output=variable,
                        sigma=model_params['sigma'], A_t=model_params['A_t'],
                        S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                        gamma_d=model_params['gamma_d'],
                        gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                        gamma_i=model_params['gamma_1'], V_d=model_params['Vd'],
                        P=model_params['P'], label=label, color=color, linestyle=linestyle,
                        linewidth=4)
    plt.ylabel('∆ energy $(pN \cdot nm)$')
    plt.legend()
    plt.tight_layout()
    plt.savefig(fig_dir+'droplet_energies_vs-H_i_delta.svg')
    plt.savefig(fig_dir+'droplet_energies_vs-H_i_delta.png')
    
  • each individual energy term output
    • set temporary variables
        P_temp = 1e-3
        gamma_d_temp = gamma_d*10
        gamma_m_temp = gamma_m*10
        K_temp = 160
      
    • normal droplet
      show_geometry(H_i=100)
      plt.savefig(fig_dir+'normal_droplet.png')
      
      • with appropriate signs

        droplet only:

        variables = ['surface_tension', 'wetting_energy', 'droplet_energy']
        colors = ['blue', 'green', 'black']
        
        plt.figure(figsize=(6,6))
        
        for variable, color in zip(variables, colors):
            plot_simulation_delta(output=variable, label=variable, P=P_temp, gamma_d=gamma_d_temp,
                                  gamma_m=gamma_m_temp, K=K_temp, color=color)
            
        plt.legend()
        #plt.ylim(-100, 100)
        plt.ylabel('∆ energy (pN * nm)')
        #plt.title(title)
        plt.tight_layout()
        plt.savefig(fig_dir+'droplet_energy_contributions.png')
        plt.savefig(fig_dir+'droplet_energy_contributions.svg')
        

        all contributions:

        variables = ['surface_tension', 'wetting_energy', 'droplet_energy',
                     'membrane_tension', 'bending_energy', 'turgor_pressure', 'F']
        colors = ['blue', 'green', 'black',
                  'red', 'purple', 'brown', 'orange']
        
        plt.figure(figsize=(6,6))
        
        for variable, color in zip(variables, colors):
            plot_simulation_delta(output=variable, label=variable,
                            sigma=model_params['sigma'], A_t=model_params['A_t'],
                            S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                            gamma_d=model_params['gamma_d'],
                            gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                            gamma_i=model_params['gamma_1'], V_d=model_params['Vd'],
                            P=model_params['P'], color=color)
        
        plt.legend()
        #plt.ylim(-100, 100)
        plt.ylabel('∆ energy (pN * nm)')
        #plt.title(title)
        plt.tight_layout()
        plt.savefig(fig_dir+'all_energy_contributions.png')
        plt.savefig(fig_dir+'all_energy_contributions.svg')
        
        variables = ['membrane_tension', 'bending_energy',
                     'surface_tension', 'wetting_energy',
                     'turgor_pressure']
        labels = ['membrane tension', 'bending energy',
                     'droplet surface tension', 'wetting energy',
                     'turgor pressure']
        colors = ['red', 'purple', 'blue', 'green', 'brown']
        
        plt.figure(figsize=(6,6))
        
        for variable, label, color in zip(variables, labels, colors):
            plot_simulation_delta(output=variable, label=label, P=1.1e-3, gamma_d=0.06,
                                  gamma_m=0.02, K=50, color=color)
        
        plt.legend()
        #plt.ylim(-100, 100)
        plt.ylabel('∆ free energy') # (pN * nm)
        #plt.title(title)
        plt.tight_layout()
        plt.savefig(fig_dir+'individual_energy_contributions.png')
        plt.savefig(fig_dir+'individual_energy_contributions.svg')
        

        effect of droplet:

        variables = ['non_droplet_energy', 'droplet_energy', 'F']
        labels = ['non-condensate energy contributions',
                  'condensate energy contributions',
                  'total energy']
        colors = ['blue', 'green', 'black']
        
        plt.figure(figsize=(6,5))
        
        for variable, label, color in zip(variables, labels, colors):
            plot_simulation_delta(output=variable, label=label, P=1.1e-3, gamma_d=0.06,
                                  gamma_m=0.02, K=50, color=color)
        
        plt.legend()
        #plt.ylim(-100, 100)
        plt.ylabel('∆ free energy ($pN \cdot nm$)') # (pN * nm)
        #plt.title(title)
        plt.tight_layout()
        plt.savefig(fig_dir+'condensate_energy_contribution.png')
        plt.savefig(fig_dir+'condensate_energy_contribution.svg')
        
      • flipped sign of negative energy contributions
        variables = ['membrane_tension', 'bending_energy_flipped','surface_tension',
                      'wetting_energy_flipped', 'turgor_pressure']
        
        plt.figure(figsize=(10,10))
        
        for variable in variables:
            plot_simulation_delta(output=variable, label=variable, P=P_temp, gamma_d=gamma_d_temp,
                                  gamma_m=gamma_m_temp, K=K_temp)
        plt.legend()
        #plt.ylim(-100, 100)
        plt.ylabel('∆ energy (pN * nm)')
        #plt.title(title)
        plt.tight_layout()
        plt.savefig(fig_dir+'energy_contributions_flipped.png')
        
    • huge droplet
      show_geometry(H_i=100, V_d=V_d*10**4,
                    left=-3600, right=3600, bottom=-2300)
      plt.savefig(fig_dir+'bigdrop.png')
      
      • with appropriate signs
        variables = ['surface_tension', 'wetting_energy', 'droplet_energy']
        colors = ['blue', 'green', 'black']
        
        plt.figure(figsize=(6,6))
        
        for variable, color in zip(variables, colors):
            plot_simulation_delta(output=variable, label=variable, P=P_temp, gamma_d=gamma_d_temp,
                                  gamma_m=gamma_m_temp, K=K_temp, V_d=V_d*10**4, color=color)
        plt.legend()
        #plt.ylim(-100, 100)
        plt.ylabel('∆ energy (pN * nm)')
        #plt.title(title)
        plt.tight_layout()
        plt.savefig(fig_dir+'droplet_energy_contributions_bigdrop.png')
        plt.savefig(fig_dir+'droplet_energy_contributions_bigdrop.svg')
        
        variables = ['non_droplet_energy', 'droplet_energy', 'F']
        labels = ['non-condensate energy contributions',
                  'condensate energy contributions',
                  'total energy']
        colors = ['blue', 'green', 'black']
        
        plt.figure(figsize=(6,5))
        
        for variable, label, color in zip(variables, labels, colors):
            plot_simulation_delta(output=variable, label=label, P=1.1e-3, gamma_d=0.06,
                                  gamma_m=0.02, K=50, V_d=V_d*10, color=color)
        
        plt.legend()
        #plt.ylim(-100, 100)
        plt.ylabel('∆ free energy ($pN \cdot nm$)') #pN*nm
        #plt.title(title)
        plt.tight_layout()
        plt.savefig(fig_dir+'condensate_energy_contribution_bigv.png')
        plt.savefig(fig_dir+'condensate_energy_contribution_bigv.svg')
        
      • flipped sign of negative energy contributions
        variables = ['membrane_tension', 'bending_energy_flipped','surface_tension',
                      'wetting_energy_flipped', 'turgor_pressure']
        
        plt.figure(figsize=(10,10))
        
        for variable in variables:
            plot_simulation_delta(output=variable, label=variable, P=P_temp, gamma_d=gamma_d_temp,
                                  gamma_m=gamma_m_temp, K=K_temp, V_d=V_d*10**4)
        plt.legend()
        #plt.ylim(-100, 100)
        plt.ylabel('∆ energy (pN * nm)')
        #plt.title(title)
        plt.tight_layout()
        plt.savefig(fig_dir+'energy_contributions_bigdrop_flipped.png')
        
  • Plot F for a series of values of gamma_m

    What are the parameter values? I assume that gamma_m is negative. Plot F for a series of values of gamma_m. Based on the plots for S_m, S_d, and S_0, F should end up with a minimum at H/R=1 if gamma_m is large enough. It might also happen if gamma_d is small enough. It’s hard to say because I can’t easily compare the magnitude of changes in S_d, S_i, and S_0.

    • absolute
           gamma_m_range = [0,50,100,150,199]
      
           for i in gamma_m_range:
               plot_simulation(gamma_m=i, label=i/gamma_d)
      
           plt.legend(title='gamma_m/gamma_d')
           plt.xlabel('invagination distance/vesicle diameter')
           plt.tight_layout()
           plt.savefig(fig_dir+'f_scan_gamma_m.png')
      
    • normalized
           gamma_m_range = [0,50,100,150,199]
      
           for i in gamma_m_range:
               plot_simulation_normalized(gamma_m=i, label=i/gamma_d)
      
           plt.legend(title='gamma_m/gamma_d')
           plt.xlabel('invagination distance/vesicle diameter')
           plt.tight_layout()
           plt.savefig(fig_dir+'f_scan_gamma_m.png')
      

      :end:

    • delta
         gamma_m_range = np.array( [0.01,0.25,0.5,0.75] )*gamma_d
      
         for i in gamma_m_range:
             plot_simulation_delta(gamma_m=i, label=i/gamma_d)
      
         plt.legend(title='$\gamma_{m}/\gamma_{d}$')
         plt.tight_layout()
         plt.savefig(fig_dir+'f_scan_gamma_m_delta.png')
      
    • delta and geometry side by side
           gamma_m_range = [0,50,100,150,180]
      
           for i in gamma_m_range:
               plot_simulation_geometry(gamma_m=i, label=i/gamma_d)
               plt.legend(title='$\gamma_{m}/\gamma_{d}$')
               plt.tight_layout()
               plt.savefig(fig_dir+'f_scan_gamma_m_%s.png' %(i))
      
  • F over series of turgor pressure
         pressure_range = np.linspace(10**-2,1,5)
    
         for i in pressure_range:
             plot_simulation_delta(P=i, label=i)
    
         plt.legend(title='turgor pressure ($N/nm^{2}$)')
         plt.tight_layout()
         plt.savefig(fig_dir+'f_scan_turgor_delta.png')
    
  • effect of droplet volume on energy barrier
    #volume_range = np.linspace(2/3*pi*200**3, 2/3*pi*2000**3,1000)
    volume_range = np.logspace(7, 11, 100)
    fc_prime_0 = Fprime(h=0, sigma=0, A_t=A_t, S_i=S_i, K=0,
               C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i,
               V_d=volume_range, P=0)
    
    plt.figure(figsize=(7,6))
    plt.plot(volume_range,fc_prime_0)
    plt.xscale('log')
    #plt.yscale('log')
    plt.xlabel('condensate volume ($nm^3$)')
    plt.ylabel(r'$\frac{dF_c}{dh} (h=0)$'+'\n (slope of condensate free energy at h=0)')
    plt.tight_layout()
    plt.savefig(fig_dir+'dfc_dh_0_vdscan.svg')
    plt.savefig(fig_dir+'dfc_dh_0_vdscan.png')
    
    #volume_range = np.linspace(2/3*pi*200**3, 2/3*pi*2000**3,1000)
    volume_range = np.logspace(7, 11, 100)
    full_fc_prime_0 = Fprime(h=0, sigma=0, A_t=A_t, S_i=S_i, K=0,
               C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i,
               V_d=volume_range, P=0)
    tension_fc_prime_0 = Fprime(h=0, sigma=0, A_t=A_t, S_i=S_i, K=0,
               C_i=C_i, gamma_d=gamma_d, gamma_m=10e-20, gamma_i=10e-20,
               V_d=volume_range, P=0)
    wetting_fc_prime_0 = Fprime(h=0, sigma=0, A_t=A_t, S_i=S_i, K=0,
               C_i=C_i, gamma_d=10e-20, gamma_m=gamma_m, gamma_i=10e-20,
               V_d=volume_range, P=0)
    
    plt.figure(figsize=(7,6))
    plt.plot(volume_range,tension_fc_prime_0,label='surface tension')
    plt.plot(volume_range,wetting_fc_prime_0,label='wetting energy')
    plt.plot(volume_range,full_fc_prime_0,label='sum')
    plt.xscale('log')
    #plt.yscale('log')
    plt.legend()
    plt.xlabel('condensate volume ($nm^3$)')
    plt.ylabel(r'$\frac{dF_c}{dh} (h=0)$'+'\n (slope of condensate free energy at h=0)')
    plt.tight_layout()
    plt.savefig(fig_dir+'dfcs_dh_0_vdscan.svg')
    plt.savefig(fig_dir+'dfcs_dh_0_vdscan.png')
    
    volume_range = np.logspace(7, 12, 100)
    dfmaxlist = []
    
    for v in volume_range:
        f = llps_model_recbrt(h=h_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                    C_i=C_i, gamma_d=0.5, gamma_m=0.2, gamma_i=0.02,
                    V_d=v, P=0.005)
        dfmax = f.max()-f[0]
        dfmaxlist.append(dfmax)
    
    plt.figure(figsize=(7,6))
    plt.plot(volume_range,dfmaxlist)
    plt.xscale('log')
    #plt.yscale('log')
    plt.xlabel('condensate volume ($nm^3$)')
    plt.ylabel(r'$\Delta F_{max}$')
    plt.tight_layout()
    plt.savefig(fig_dir+'dfmax_vdscan.svg')
    plt.savefig(fig_dir+'dfmax_vdscan.png')
    
    • plot energy derivatives vs. volume

      droplet only

      V_d_range = np.linspace(5E6,5E7)
      
      plt.plot(V_d_range,Fprime(h=0,
                          sigma=0, A_t=model_params['A_t'],
                          S_i=model_params['A_c'], K=0, C_i=model_params['H1'],
                          gamma_d=model_params['gamma_d'],
                          gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                          gamma_i=model_params['gamma_1'], V_d=V_d_range,
                                P=0), label=r'$\frac{d F_{d}}{d h}$', color='black',# linestyle=linestyle,
                          linewidth=4)
      
      # plt.plot(V_d_range,Fprime(h=0,
      #                     sigma=model_params['sigma'], A_t=model_params['A_t'],
      #                     S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
      #                     gamma_d=model_params['gamma_d'],
      #                     gamma_m=model_params['gamma_3']-model_params['gamma_2'],
      #                     gamma_i=model_params['gamma_1'], V_d=V_d_range,
      #                     P=model_params['P']), label='F all', color='orange', linestyle=linestyle,
      #                     linewidth=4)
      
      plt.plot(V_d_range,-Fprime(h=0,
                          sigma=model_params['sigma'], A_t=model_params['A_t'],
                          S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                          gamma_d=model_params['gamma_d'],
                          gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                          gamma_i=model_params['gamma_1'], V_d=V_d_range,
                          P=model_params['P']) +
               Fprime(h=0,
                          sigma=0, A_t=model_params['A_t'],
                          S_i=model_params['A_c'], K=0, C_i=model_params['H1'],
                          gamma_d=model_params['gamma_d'],
                          gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                          gamma_i=model_params['gamma_1'], V_d=V_d_range,
                                P=0), label=r'$-\frac{d F_{e}}{d h}$', color='orange',# linestyle=linestyle,
                          linewidth=4)
      
      plt.legend()
      plt.xlabel('$V_d$: droplet volume')
      plt.ylabel('slope of free energy at $h=0$')
      plt.tight_layout()
      plt.savefig(fig_dir+'dfdh_vdscan.svg')
      plt.savefig(fig_dir+'dfdh_vdscan.png')
      
    • animation
      animate_simulation_geometry(
                          sigma=model_params['sigma'], A_t=model_params['A_t']*1000,
                          S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                          gamma_d=model_params['gamma_d'],
                          gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                          gamma_i=model_params['gamma_1'], V_d_range=V_d_range,
                          P=model_params['P'],
                      geometry_top=20, geometry_bottom=-270, geometry_left=-400, geometry_right=400,
                      plot_top=2300, plot_bottom=-1900, plot_left=-0.1, plot_right=0.1
                          )
      
  • spontaneity analysis

    there’s no solution for kinetic \(\sigma\)

    • DONE \(V_{d}\) x-axis
      boundaries = [(thermo_spont_K,kinetic_spont_K, 'K', 2),
                    (thermo_spont_P,kinetic_spont_P, 'P', -2),
                    #(thermo_spont_sigma,kinetic_spont_sigma, 'sigma'),
                    (thermo_spont_C_i,kinetic_spont_C_i, 'C_i', 2)]
      x_range = np.logspace(7,9,base=10)
      for thermo_boundaryy, kinetic_boundaryy, y_param, y_offset in boundaries:
          thermo_boundary = thermo_boundaryy(V_d=x_range)
          kinetic_boundary = kinetic_boundaryy(V_d=x_range)
          both_boundary_test(x_param='V_d', y_param=y_param, x_range=x_range,
                          thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                          y_offset=y_offset/2, factors=[0.9,1,1.1])
          plt.savefig(fig_dir+y_param+'V_d_boundary_check.png')
      
      show_geometry(H_i=100, V_d=10**(7))
      plt.savefig(fig_dir+'10e7_geometry.png')
      
      show_geometry(H_i=100, V_d=10**(9), left=-1200, right=1200, bottom=-600)
      plt.savefig(fig_dir+'10e9_geometry.png')
      
    • DONE \(\gamma_{d}\) x-axis
      boundaries = [(thermo_spont_K,kinetic_spont_K, 'K', -4),
                    (thermo_spont_P,kinetic_spont_P, 'P', 2),
                    #(thermo_spont_sigma,kinetic_spont_sigma, 'sigma'),
                    (thermo_spont_C_i,kinetic_spont_C_i, 'C_i', 2)]
      x_range = np.logspace(-1,0,base=10)
      for thermo_boundaryy, kinetic_boundaryy, y_param, y_offset in boundaries:
          thermo_boundary = thermo_boundaryy(gamma_d=x_range)
          kinetic_boundary = kinetic_boundaryy(gamma_d=x_range)
          both_boundary_test(x_param='gamma_d', y_param=y_param, x_range=x_range,
                          thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                          y_offset=y_offset/2, factors=[0.95,1,1.05])
          plt.savefig(fig_dir+y_param+'gamma_d_boundary_check.png')
          plt.savefig(fig_dir+y_param+'gamma_d_boundary_check.svg')
      
    • DONE \(\gamma_{m}\) x-axis
      boundaries = [(thermo_spont_K,kinetic_spont_K, 'K', 2),
                    (thermo_spont_P,kinetic_spont_P, 'P', -2),
                    #(thermo_spont_sigma,kinetic_spont_sigma, 'sigma'),
                    (thermo_spont_C_i,kinetic_spont_C_i, 'C_i', 2)]
      x_range = np.logspace(-3,0,base=10)
      for thermo_boundaryy, kinetic_boundaryy, y_param, y_offset in boundaries:
          thermo_boundary = thermo_boundaryy(gamma_m=x_range)
          kinetic_boundary = kinetic_boundaryy(gamma_m=x_range)
          both_boundary_test(x_param='gamma_m', y_param=y_param, x_range=x_range,
                          thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                          y_offset=y_offset*2, factors=[0.95,1,1.05])
          plt.savefig(fig_dir+y_param+'gamma_m_boundary_check.png')
      
    • TODO \(\gamma_{m}\) /\(\gamma_{d}\) x-axis
      boundaries = [(thermo_spont_K,kinetic_spont_K, 'K', 2),
                    (thermo_spont_P,kinetic_spont_P, 'P', -2),
                    #(thermo_spont_sigma,kinetic_spont_sigma, 'sigma'),
                    (thermo_spont_C_i,kinetic_spont_C_i, 'C_i', 2)]
      x_range = np.logspace(-2,0.0001,base=10)*gamma_d
      for thermo_boundaryy, kinetic_boundaryy, y_param, y_offset in boundaries:
          thermo_boundary = thermo_boundaryy(gamma_m=x_range)
          kinetic_boundary = kinetic_boundaryy(gamma_m=x_range)
          both_boundary_test(x_param='gamma_m', y_param=y_param, x_range=x_range,
                          thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                          y_offset=y_offset/2, factors=[0.8,1,1.2])
          plt.savefig(fig_dir+y_param+'gamma_m_relative_boundary_check.png')
      
    • DONE \(\sigma\) x-axis
      boundaries = [(thermo_spont_K,kinetic_spont_K, 'K', -2),
                    (thermo_spont_P,kinetic_spont_P, 'P', 2),
                    (thermo_spont_C_i,kinetic_spont_C_i, 'C_i', -2)]
      x_range = np.logspace(-5,-1,base=10)
      for thermo_boundaryy, kinetic_boundaryy, y_param, y_offset in boundaries:
          thermo_boundary = thermo_boundaryy(sigma=x_range)
          kinetic_boundary = kinetic_boundaryy(sigma=x_range)
          both_boundary_test(x_param='sigma', y_param=y_param, x_range=x_range,
                          thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                          y_offset=y_offset/5, factors=[0.8,1,1.2])
          plt.savefig(fig_dir+y_param+'sigma_boundary_check.png')
      
    • DONE \(S_{i}\) x-axis
      boundaries = [(thermo_spont_K,kinetic_spont_K, 'K', 2),
                    (thermo_spont_P,kinetic_spont_P, 'P', -2),
                    (thermo_spont_C_i,kinetic_spont_C_i, 'C_i', 2)]
      x_range = np.logspace(3.5,4.5,base=10)
      for thermo_boundaryy, kinetic_boundaryy, y_param, y_offset in boundaries:
          thermo_boundary = thermo_boundaryy(S_i=x_range)
          kinetic_boundary = kinetic_boundaryy(S_i=x_range)
          both_boundary_test(x_param='S_i', y_param=y_param, x_range=x_range,
                          thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                          y_offset=y_offset*50, factors=[0.999,1,1.001])
          plt.savefig(fig_dir+y_param+'S_i_boundary_check.png')
      
    • DONE \(P\) x-axis
      boundaries = [(thermo_spont_K,kinetic_spont_K, 'K', 2),
                    #(thermo_spont_P,kinetic_spont_P, 'P', 2),
                    (thermo_spont_C_i,kinetic_spont_C_i, 'C_i', 2)]
      x_range = np.logspace(-3,0,base=10)
      for thermo_boundaryy, kinetic_boundaryy, y_param, y_offset in boundaries:
          thermo_boundary = thermo_boundaryy(P=x_range)
          kinetic_boundary = kinetic_boundaryy(P=x_range)
          both_boundary_test(x_param='P', y_param=y_param, x_range=x_range,
                          thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                          y_offset=y_offset*150, factors=[0.999,1,1.001])
          plt.savefig(fig_dir+y_param+'P_boundary_check.png')
      
    • DONE \(C_{i}\) x-axis
      boundaries = [(thermo_spont_K,kinetic_spont_K, 'K', 2),
                    (thermo_spont_P,kinetic_spont_P, 'P', -2)]
      x_range = np.logspace(-4,-1,base=10)
      for thermo_boundaryy, kinetic_boundaryy, y_param, y_offset in boundaries:
          thermo_boundary = thermo_boundaryy(C_i=x_range)
          kinetic_boundary = kinetic_boundaryy(C_i=x_range)
          both_boundary_test(x_param='C_i', y_param=y_param, x_range=x_range,
                          thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                          y_offset=y_offset/5, factors=[0.8,1,1.2])
          plt.savefig(fig_dir+y_param+'C_i_boundary_check.png')
      
    • DONE \(K\) x-axis
      boundaries = [(thermo_spont_C_i,kinetic_spont_C_i, 'C_i', 2),
                    (thermo_spont_P,kinetic_spont_P, 'P', -2)]
      x_range = np.logspace(1,4,base=10)
      for thermo_boundaryy, kinetic_boundaryy, y_param, y_offset in boundaries:
          thermo_boundary = thermo_boundaryy(K=x_range)
          kinetic_boundary = kinetic_boundaryy(K=x_range)
          both_boundary_test(x_param='K', y_param=y_param, x_range=x_range,
                          thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                          y_offset=y_offset/1, factors=[0.9,1,1.1])
          plt.savefig(fig_dir+y_param+'K_boundary_check.png')
      
  • plots for 20200915 lab meeting drafting
    • TODO each individual energy term output
          variables = ['membrane_tension',
                        'bending_energy','surface_tension',
                        'wetting_energy','turgor_pressure']
      
          for variable in variables:
              plt.figure()
              plot_simulation_delta(output=variable)
              plt.tight_layout()
              plt.savefig(fig_dir+'%s-vs-H_i_delta.png' %(variable))
      
    • membrane components without droplet
      • effect of membrane tension alone
         plot_simulation_geometry(V_d=0, K=0, gamma_m=0, gamma_i=0, P=0)
        
      • effect of bending energy alone
         plot_simulation_geometry(K=K, V_d=0, sigma=0, gamma_m=0, gamma_i=0, P=0)
        
        • TODO fix this

          why is this different from plotting membrane bending energy term alone?

           R_i_range = np.linspace(0.1,np.pi,100)
           bending_energy = (K/2)*S_i*((2/R_i_range)-C_i)**2
           plt.plot(R_i_range, bending_energy)
          
    • droplet effects
      • surface tension
         plot_simulation_geometry(sigma=0, K=0, gamma_m=0, gamma_i=0, P=0)
        
      • surface tension + wetting energy

        can’t have wetting energy alone

         plot_simulation_geometry(sigma=0, K=0, P=0)
        
    • TODO droplet+membrane effects
    • TODO turgor pressure effects
  • complementing membrane mechanics model
    • kinetically spontaneous
      N=10
      #coat_area = np.pi*100**2
      #coat_radius = np.sqrt(coat_area/(4*np.pi))
      #total_radius = np.sqrt(coat_area/np.pi)*3
      #Vd = 1/3*np.pi*200**3
      
      coat_radius = 50
      #coat_area = ((3*N)/(3*N+2))*4*np.pi*coat_radius**2
      coat_area = 4*np.pi*coat_radius**2
      total_radius = np.sqrt(coat_area/np.pi)*3
      Vd = 1/2*np.pi*200**3
      
      model_params = {
          'time':1000000,
          'dt':1,
          'print_interval':10000,
          'write_interval':1000,
          'shape':'initiation',
          'P':0.0001,
          'sigma':0.01,
          'gamma_1':0,
          'gamma_2':0.002,
          'gamma_3':0.02,
          'gamma_d':0.04,
          'h1':1E-2,
          'h2':1E-2,
          'h3':1E-2,
          'H1':1/coat_radius,
          'K1':50,
          'H3':0,
          'K3':1,
          'polar_type':'integral',
          'G':1,
          'cell_wall':True,
          'N':N,
          'A_t':total_radius,
          'A_c':coat_area,
          'Vd':Vd,
          'emap':False,
          'tip_pull':0
      }
      
      
      plot_simulation_geometry(x=H_i_range/vesicle_diameter, output='F',
                           H_i=H_i_range, sigma=model_params['sigma'], A_t=model_params['A_t'],
                          S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                          gamma_d=model_params['gamma_d'], gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                          gamma_i=model_params['gamma_1'], V_d=model_params['Vd'], P=model_params['P'],
                           fig_width=15, fig_height=5, label=None, variable_df=variable_df,
                          num_plots=2, top=20, bottom=-250, left=-400, right=400)
      
      variables = ['surface_tension', 'wetting_energy', 'droplet_energy',
                   'membrane_tension', 'bending_energy', 'turgor_pressure', 'F']
      colors = ['blue', 'green', 'black',
                'red', 'purple', 'brown', 'orange']
      
      plt.figure(figsize=(6,6))
      
      for variable, color in zip(variables, colors):
          plot_simulation_delta(output=variable, label=variable, color=color,
                           H_i=H_i_range, sigma=model_params['sigma'], A_t=model_params['A_t'],
                          S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                          gamma_d=model_params['gamma_d'], gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                          gamma_i=model_params['gamma_1'], V_d=model_params['Vd'], P=model_params['P'],
                           variable_df=variable_df)
      
      plt.legend()
      plt.hlines(0,-0.1,1.1,linestyles='dotted',colors='black')
      #plt.ylim(-100, 100)
      plt.ylabel('∆ energy (pN * nm)')
      #plt.title(title)
      plt.tight_layout()
      plt.savefig(fig_dir+'energy_contributions_complementary.png')
      plt.savefig(fig_dir+'energy_contributions_complementary.svg')
      
      plot_simulation_delta(output='F', label='F', color='orange',
                           H_i=H_i_range, sigma=model_params['sigma'], A_t=model_params['A_t'],
                          S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                          gamma_d=model_params['gamma_d'], gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                          gamma_i=model_params['gamma_1'], V_d=model_params['Vd'], P=model_params['P'],
                           variable_df=variable_df)
      plt.legend()
      #plt.hlines(0,-0.1,1.1,linestyles='dotted',colors='black')
      #plt.ylim(-100, 100)
      plt.ylabel('∆ energy (pN * nm)')
      #plt.title(title)
      plt.tight_layout()
      plt.savefig(fig_dir+'f_complementary.png')
      plt.savefig(fig_dir+'f_complementary.svg')
      
    • only thermodynamically spontaneous
      N=10
      #coat_area = np.pi*100**2
      #coat_radius = np.sqrt(coat_area/(4*np.pi))
      #total_radius = np.sqrt(coat_area/np.pi)*3
      #Vd = 1/3*np.pi*200**3
      
      coat_radius = 50
      #coat_area = ((3*N)/(3*N+2))*4*np.pi*coat_radius**2
      coat_area = 4*np.pi*coat_radius**2
      total_radius = np.sqrt(coat_area/np.pi)*3
      Vd = 1/2*np.pi*200**3
      
      model_params = {
          'time':1000000,
          'dt':1,
          'print_interval':10000,
          'write_interval':1000,
          'shape':'initiation',
          'P':0.0001,
          'sigma':0.01,
          'gamma_1':0,
          'gamma_2':0.002,
          'gamma_3':0.02,
          'gamma_d':0.06,
          'h1':1E-2,
          'h2':1E-2,
          'h3':1E-2,
          'H1':1/coat_radius,
          'K1':10,
          'H3':0,
          'K3':1,
          'polar_type':'integral',
          'G':1,
          'cell_wall':True,
          'N':N,
          'A_t':total_radius,
          'A_c':coat_area,
          'Vd':Vd,
          'emap':False,
          'tip_pull':0
      }
      
      
      plot_simulation_geometry(x=H_i_range/vesicle_diameter, output='F',
                           H_i=H_i_range, sigma=model_params['sigma'], A_t=model_params['A_t'],
                          S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                          gamma_d=model_params['gamma_d'], gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                          gamma_i=model_params['gamma_1'], V_d=model_params['Vd'], P=model_params['P'],
                           fig_width=15, fig_height=5, label=None, variable_df=variable_df,
                          num_plots=2, top=20, bottom=-250, left=-400, right=400)
      
      variables = ['surface_tension', 'wetting_energy', 'droplet_energy',
                   'membrane_tension', 'bending_energy', 'turgor_pressure', 'F']
      colors = ['blue', 'green', 'black',
                'red', 'purple', 'brown', 'orange']
      
      plt.figure(figsize=(6,6))
      
      for variable, color in zip(variables, colors):
          plot_simulation_delta(output=variable, label=variable, color=color,
                           H_i=H_i_range, sigma=model_params['sigma'], A_t=model_params['A_t'],
                          S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                          gamma_d=model_params['gamma_d'], gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                          gamma_i=model_params['gamma_1'], V_d=model_params['Vd'], P=model_params['P'],
                           variable_df=variable_df)
      
      plt.legend()
      plt.hlines(0,-0.1,1.1,linestyles='dotted',colors='black')
      #plt.ylim(-100, 100)
      plt.ylabel('∆ energy (pN * nm)')
      #plt.title(title)
      plt.tight_layout()
      plt.savefig(fig_dir+'energy_contributions_complementary_barrier.png')
      plt.savefig(fig_dir+'energy_contributions_complementary_barrier.svg')
      
      plot_simulation_delta(output='F', label='F', color='orange',
                           H_i=H_i_range, sigma=model_params['sigma'], A_t=model_params['A_t'],
                          S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                          gamma_d=model_params['gamma_d'], gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                          gamma_i=model_params['gamma_1'], V_d=model_params['Vd'], P=model_params['P'],
                           variable_df=variable_df)
      plt.legend()
      #plt.hlines(0,-0.1,1.1,linestyles='dotted',colors='black')
      #plt.ylim(-100, 100)
      plt.ylabel('∆ energy (pN * nm)')
      #plt.title(title)
      plt.tight_layout()
      plt.savefig(fig_dir+'f_complementary_barrier.png')
      plt.savefig(fig_dir+'f_complementary_barrier.svg')
      
    • not spontaneous
      N=10
      #coat_area = np.pi*100**2
      #coat_radius = np.sqrt(coat_area/(4*np.pi))
      #total_radius = np.sqrt(coat_area/np.pi)*3
      #Vd = 1/3*np.pi*200**3
      
      coat_radius = 50
      #coat_area = ((3*N)/(3*N+2))*4*np.pi*coat_radius**2
      coat_area = 4*np.pi*coat_radius**2
      total_radius = np.sqrt(coat_area/np.pi)*3
      Vd = 1/2*np.pi*200**3
      
      model_params = {
          'time':1000000,
          'dt':1,
          'print_interval':10000,
          'write_interval':1000,
          'shape':'initiation',
          'P':0.0001,
          'sigma':0.01,
          'gamma_1':0,
          'gamma_2':0.002,
          'gamma_3':0.02,
          'gamma_d':0.06,
          'h1':1E-2,
          'h2':1E-2,
          'h3':1E-2,
          'H1':1/coat_radius,
          'K1':1,
          'H3':0,
          'K3':1,
          'polar_type':'integral',
          'G':1,
          'cell_wall':True,
          'N':N,
          'A_t':total_radius,
          'A_c':coat_area,
          'Vd':Vd,
          'emap':False,
          'tip_pull':0
      }
      
      
      plot_simulation_geometry(x=H_i_range/vesicle_diameter, output='F',
                           H_i=H_i_range, sigma=model_params['sigma'], A_t=model_params['A_t'],
                          S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                          gamma_d=model_params['gamma_d'], gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                          gamma_i=model_params['gamma_1'], V_d=model_params['Vd'], P=model_params['P'],
                           fig_width=15, fig_height=5, label=None, variable_df=variable_df,
                          num_plots=2, top=20, bottom=-250, left=-400, right=400)
      
      variables = ['surface_tension', 'wetting_energy', 'droplet_energy',
                   'membrane_tension', 'bending_energy', 'turgor_pressure', 'F']
      colors = ['blue', 'green', 'black',
                'red', 'purple', 'brown', 'orange']
      
      plt.figure(figsize=(6,6))
      
      for variable, color in zip(variables, colors):
          plot_simulation_delta(output=variable, label=variable, color=color,
                           H_i=H_i_range, sigma=model_params['sigma'], A_t=model_params['A_t'],
                          S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                          gamma_d=model_params['gamma_d'], gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                          gamma_i=model_params['gamma_1'], V_d=model_params['Vd'], P=model_params['P'],
                           variable_df=variable_df)
      
      plt.legend()
      plt.hlines(0,-0.1,1.1,linestyles='dotted',colors='black')
      #plt.ylim(-100, 100)
      plt.ylabel('∆ energy (pN * nm)')
      #plt.title(title)
      plt.tight_layout()
      plt.savefig(fig_dir+'energy_contributions_complementary_notspon.png')
      plt.savefig(fig_dir+'energy_contributions_complementary_notspon.svg')
      
      plot_simulation_delta(output='F', label='F', color='orange',
                           H_i=H_i_range, sigma=model_params['sigma'], A_t=model_params['A_t'],
                          S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                          gamma_d=model_params['gamma_d'], gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                          gamma_i=model_params['gamma_1'], V_d=model_params['Vd'], P=model_params['P'],
                           variable_df=variable_df)
      plt.legend()
      #plt.hlines(0,-0.1,1.1,linestyles='dotted',colors='black')
      #plt.ylim(-100, 100)
      plt.ylabel('∆ energy (pN * nm)')
      #plt.title(title)
      plt.tight_layout()
      plt.savefig(fig_dir+'f_complementary_notspon.png')
      plt.savefig(fig_dir+'f_complementary_notspon.svg')
      
    • kinetically spontaneous flat
      N=10
      #coat_area = np.pi*100**2
      #coat_radius = np.sqrt(coat_area/(4*np.pi))
      #total_radius = np.sqrt(coat_area/np.pi)*3
      #Vd = 1/3*np.pi*200**3
      
      coat_radius = 50
      #coat_area = ((3*N)/(3*N+2))*4*np.pi*coat_radius**2
      coat_area = 4*np.pi*coat_radius**2
      total_radius = np.sqrt(coat_area/np.pi)*3
      Vd = 1/2*np.pi*200**3
      
      model_params = {
          'time':1000000,
          'dt':1,
          'print_interval':10000,
          'write_interval':1000,
          'shape':'initiation',
          'P':0.0001,
          'sigma':0.1,
          'gamma_1':0,
          'gamma_2':0.002,
          'gamma_3':0.02,
          'gamma_d':0.06,
          'h1':1E-2,
          'h2':1E-2,
          'h3':1E-2,
          'H1':1/coat_radius,
          'K1':1,
          'H3':0,
          'K3':1,
          'polar_type':'integral',
          'G':1,
          'cell_wall':True,
          'N':N,
          'A_t':total_radius,
          'A_c':coat_area,
          'Vd':Vd,
          'emap':False,
          'tip_pull':0
      }
      
      
      plot_simulation_geometry(x=H_i_range/vesicle_diameter, output='F',
                           H_i=H_i_range, sigma=model_params['sigma'], A_t=model_params['A_t'],
                          S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                          gamma_d=model_params['gamma_d'], gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                          gamma_i=model_params['gamma_1'], V_d=model_params['Vd'], P=model_params['P'],
                           fig_width=15, fig_height=5, label=None, variable_df=variable_df,
                          num_plots=2, top=20, bottom=-250, left=-400, right=400)
      
      variables = ['surface_tension', 'wetting_energy', 'droplet_energy',
                   'membrane_tension', 'bending_energy', 'turgor_pressure', 'F']
      colors = ['blue', 'green', 'black',
                'red', 'purple', 'brown', 'orange']
      
      plt.figure(figsize=(6,6))
      
      for variable, color in zip(variables, colors):
          plot_simulation_delta(output=variable, label=variable, color=color,
                           H_i=H_i_range, sigma=model_params['sigma'], A_t=model_params['A_t'],
                          S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                          gamma_d=model_params['gamma_d'], gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                          gamma_i=model_params['gamma_1'], V_d=model_params['Vd'], P=model_params['P'],
                           variable_df=variable_df)
      
      plt.legend()
      plt.hlines(0,-0.1,1.1,linestyles='dotted',colors='black')
      #plt.ylim(-100, 100)
      plt.ylabel('∆ energy (pN * nm)')
      #plt.title(title)
      plt.tight_layout()
      plt.savefig(fig_dir+'energy_contributions_complementary_sponflat.png')
      plt.savefig(fig_dir+'energy_contributions_complementary_sponflat.svg')
      
      plot_simulation_delta(output='F', label='F', color='orange',
                           H_i=H_i_range, sigma=model_params['sigma'], A_t=model_params['A_t'],
                          S_i=model_params['A_c'], K=model_params['K1'], C_i=model_params['H1'],
                          gamma_d=model_params['gamma_d'], gamma_m=model_params['gamma_3']-model_params['gamma_2'],
                          gamma_i=model_params['gamma_1'], V_d=model_params['Vd'], P=model_params['P'],
                           variable_df=variable_df)
      plt.legend()
      #plt.hlines(0,-0.1,1.1,linestyles='dotted',colors='black')
      #plt.ylim(-100, 100)
      plt.ylabel('∆ energy (pN * nm)')
      #plt.title(title)
      plt.tight_layout()
      plt.savefig(fig_dir+'f_complementary_sponflat.png')
      plt.savefig(fig_dir+'f_complementary_sponflat.svg')
      

Author: Max Ferrin

Created: 2024-04-02 Tue 11:00